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The  accomplishments  of  the  research  project  are  advancements  in  the  state  of 
knowledge  in  the  reliability  assessment  of  ship  structures  against  fatigue. 
Fatigue  life  of  ship  structural  components  is  believed  sensitive  to  load  events 
exceeding  certain  thresholds.  Ocean  waves  are  known  to  be  the  dominant 
source  of  fatigue  in  ship  structures.  Thus,  this  fatigue  investigation  has  spanned 
three  technical  areas:  direct  prediction  of  fatigue  damage  in  non-Gaussian  load 
environments,  prediction  of  probability  of  exceeding  specified  load  thresholds  in 
non-Gaussian  environments,  and  simulation  and  prediction  of  non-linear,  non- 
Gaussian  ocean  waves.  Some  of  the  methods  developed  are  applicable  to  many 
different  structural  types  beyond  ship  structural  components.  In  these  cases,  the 
methods  have  been  additionally  tested  and  verified  against  wind-turbine  fatigue 
data  because  turbine-blade  fatigue  data  was  available  to  the  project. 

In  the  predicting  fatigue  damage  in  non-Gaussian  load  environments,  research 
has  concentrated  on  the  fitting  of  up  to  three  statistical  models  to  observed  load 
data  and  on  application  of  these  models  to  predict  long-term  fatigue 
accumulation  from  short-term  observed  data.  A  broad  variety  of  statistical 
models  have  been  fit  using  the  method  of  moments,  and  the  resulting  models 
have  been  tested  against  measured  data  for  ship  fatigue  and  for  wind-turbine 
blade  fatigue  applications. 

In  the  prediction  of  exceeding  specified  load  thresholds,  use  of  four-moment 
Hermite  models,  analytic  formulae,  first-order  reliability  method  (FORM), 
statistical  bootstrapping  methods,  and  fitting  of  various  statistical  models  have  all 
been  proposed  and  tested  against  measured  physical  data.  Numerical 
algorithms  in  the  form  of  computer  software  have  been  developed  as  in  an  effort 
to  more  fully  understand  these  phenomena.  Example  applications  have  included 
prediction  of  extreme  statistics  of  wind-driven  ocean  waves,  and  of  motions  of  a 
spar  production  platform  subject  to  these  waves. 

In  the  prediction  of  non-linear  non-Gaussian  ocean  waves,  a  methodology  has 
been  developed  to  simulate  second-order  random  ocean  waves  to  match  a  target 
sea  spectrum.  The  inverse  has  also  been  developed:  identifying  the  first-  and 
second-order  components  underlying  a  specified  wave  time-history.  The  results 
of  both  the  simulation  and  identification  methodologies  can  be  used  to  predict 
consistent  wave  time-histories  at  other  spatial  locations.  Results  of  this 
methodology  have  been  verified  against  observed  wave  time-histories  both  from 
field  data  and  from  model  test  measurements. 
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Abstract 


The  main  objective  of  this  thesis  is  to  study  procedures  to  estimate  statistics  of  extreme 
responses  to  random  excitation.  Various  methods  are  applied  to  predict  the  extreme 
statistics  of  the  horizontal  offset  of  a  spar  buoy,  consisting  of  3  dominant  frequency 
bands;  resonant  surge,  resonant  pitch,  and  wave  frequency  response.  In  order  to  study  the 
effect  of  limited  data,  both  1  hour  of  model  test  data  and  36  hours  of  computer  simulated 
data  were  analyzed.  Gumbel,  standard  Weibull,  quadratic  Weibull,  and  Hermite  models 
were  fitted  using  the  method  of  moments.  The  local  peaks  of  the  components  of  the 
response  appeared  to  be  well  described  by  the  Rayleigh  model  (of  which  the  Weibull  is  a 
generalization)  when  fitted  to  the  simulated  data.  The  distribution  of  the  total  response, 
being  a  mixture  of  3  Rayleigh  distributed  variables,  was  better  described  by  the  Hermite 
(local  peaks)  and  Gumbel  distributions  (global  peaks).  For  the  observed  1  hour  of  data, 
the  distributions  of  the  components  seemed  considerably  more  narrow  and  “pinched”.  We 
show  this  to  be  largely  due  to  the  effect  of  limited  data,  which  biases  the  higher  moments 
in  particular.  However,  it  is  not  unlikely  that  limitations  of  the  computational  model 
account  partly  for  this  discrepancy.  Estimates  of  the  uncertainty  in  our  predictions  were 
made  using  non-parametric  bootstrapping,  which  were  compared  to  results  from 
simulation.  The  results  suggest  that  computer  simulations  provide  good  results,  while 
non-parametric  bootstrapping  appears  to  be  less  suitable  to  estimate  statistics  of  extremes. 
In  addition,  it  appears  that  for  all  models  the  uncertainty  in  the  predictions  of  the  mean 
max  can  be  reduced  by  dividing  the  limited  data  into  short-duration  segments,  from 
which  many  estimates  of  the  mean  max  can  be  made  by  extrapolating  a  fitted  model. 
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Introduction 


Objective 

The  main  objective  of  this  thesis  is  to  study  procedures  to  estimate  statistics  of  extreme 
responses  to  random,  excitation.  The  specific  application  here  is  to  the  response  of 
floating  structures,  though  the  methods  surveyed  here  have  more  general  applicability. 
With  the  exploration  of  oil  and  gas  in  very  deep  water,  such  as  the  Gulf  of  Mexico,  West 
Africa,  and  the  North  Sea,  floating  structures  have  become  a  very  important  concept  for 
operators,  replacing  jacket  and  tower  stmctures.  In  deep  water,  jacket  and  tower 
structures  are  at  a  disadvantage  as  they  become  too  flexible.  Their  natural  frequencies 
decrease  towards  values  where  there  is  a  significant  amount  of  wave  energy,  causing 
dynamic  excitation.  Floating  structures  generally  have  their  natural  frequencies  of 
oscillation  well  away  from  the  primary  wave  frequencies.  Although  they  can  still  be 
excited  dynamically,  they  are  considerably  more  cost  effective  in  deeper  water. 

While  wave  loads  on  a  jacket  or  tower  structure  can  be  estimated  in  a  fairly 
straightforward  manner  with  Morison’s  empirical  formula,  the  loads  on  floaters  generally 
require  more  complicated  diffraction  models.  An  important  complication  is  that  the 
dynamic  response  of  floaters  is  governed  by  two  random  variables,  the  wave  height  and 
the  wave  period,  while  in  the  case  of  a  regular  jacket  or  tower  structure  only  the  wave 
height  dominates  the  quasi-static  response.  The  use  of  more  complicated  (non-linear) 
response  models  complicates  the  statistical  analysis  of  the  response.  Various  methods, 
which  can  be  used  in  the  design  process  will  be  discussed  and  applied  for  a  spar  buoy  in 


1 


this  report.  The  accurate  prediction  of  extreme  behavior  is  critical  in  the  design  process. 
The  shift  to  floating  structures  has  created  a  need  for  appropriate  statistical  models.  This 
work  is  intended  to  provide  some  insight  and  explore  some  possibilities  to  address  this 
need. 

Background 

Limit  States 

During  the  design  of  a  structure  one  has  to  consider  all  possible  ways  the  structure  can 
fail,  and  assure  that  the  probability  of  failure  is  sufficiently  small.  Failure  implies  that  the 
structure  cannot  satisfy  its  functional  requirements  anymore,  which  need  not  be 
equivalent  to  collapse.  These  functional  requirements  are  generally  modeled  as  limit 
states  Mi(xi,X2,...,Xn),  which  can  be  functions  of  many  random  variables.  The  probability 
that  a  certain  limit  state  is  exceeded  (Mi  <  0)  determines  how  reliable  or  safe  the  structure 
is  with  regard  to  that  limit  state.  The  combination  of  all  limit  states  determines  the  total 
probability  of  failure.  Generally  a  limit  state  is  expressed  as  M  =  R  -  L,  where  R 
represents  the  resistance,  and  L  represents  the  load.  More  generally,  R  and  L  should  be 
considered  general  “capacity”  and  “demand”  variables  respectively.  These  may  be 
expressed  in  terms  of  loads,  load  effects,  or  resulting  response  quantities.  Many  reliability 
studies  now  consider  displacement  (or  ductility)  demands,  as  we  will  do  here.  For 
offshore  structures  subject  to  wave  loads,  the  analysis  generally  comprises  of  two  time- 
scales,  which  are  typically  referred  to  as  “long  term”  and  “short  term”.  The  motivation 
for  this  divide  is  that  the  wave  environment,  and  hence  the  structural  response,  is 
generally  a  slowly  evolving  non-stationary  random  process.  Generally  this  situation  is 


modeled  by  assuming  the  random  wave  process  to  be  stationary  for  short-term  periods  or 
sea  states,  and  calculating  the  statistics  for  a  sea  state  conditional  on  the  long-term 
environmental  parameters,  by  which  it  is  defined.  Both  time  scales  are  discussed  in  the 
following  paragraph. 


The  Extreme  Response 

The  limit  state  we  are  concerned  with  here,  is  failure  caused  by  extreme  overload  due  to 
ocean  waves.  The  particular  case  we  will  analyze  in  chapters  3  and  4  is  the  extreme 
horizontal  offset  of  a  spar  buoy  (where,  as  noted  above,  the  offset  is  the  load  L).  If  the 
spar  exceeds  a  certain  limiting  offset  (i.e.  the  resistance  R),  failure  occurs.  In  many  cases 
the  variability  in  the  resistance  is  negligible  compared  to  the  variability  in  the  load,  and 
the  problem  reduces  to  the  load  exceeding  a  certain  deterministic  value.  The  random 
variables  that  determine  the  statistics  of  the  waves  (and  the  response)  are  the  significant 
wave  height  Hs,  the  spectral  peak  period  Tp,  and  the  extreme  response  Xmax  given  these 
two  parameters  that  define  the  stationary  random  process.  The  quantity  Hs  determines  the 
area  Gr^  =  (Hs/4)  ,  under  the  wave  power  spectrum,  while  Tp  corresponds  to  the 
frequency  with  maximum  power.  The  so  called  “short  term”  problem  is  to  estimate  the 
statistics  of  Xma,  given  Hs  and  Tp,  while  the  “long  term”  problem  then  includes 
randomness  in  Hs  and  Tp  across  various  sea  states.  We  focus  here  on  the  “short  term” 
problem  only,  assuming  Hs  and  Tp  are  fixed  (coinciding,  for  example,  with  available 
steady  state  model  test  conditions).  The  long  term  problem  can  be  handled  by 
conventional  methods;  e.g.  Form  (Madsen,  1986)  or  simulation  (Melchers,  1987). 
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The  extreme  response,  Xmax,  will  be  a  random  variable,  which  can  be  modeled,  by 
considering  the  response  process  x(t)  over  different  time  scales.  The  various  time  scales, 
and  appropriate  distribution  models  for  the  standard  Gaussian  case  are  discussed  below. 
We  begin  here  at  the  finest  time  scale,  and  proceed  to  increasingly  global  time  scales. 

Model  of  the  Entire  Process,  X(t). 

At  the  finest  time  scale,  we  may  seek  to  model  the  cumulative  distribution  function 
(CDF)  Fx(x)  of  the  random  process  x(t)  selected  at  arbitrary  time  t: 

F,M=P[x{l)ix]  l-I 

In  the  most  common  case  X(t)  is  assumed  Gaussian,  in  which  case  Fx(x)  can  be  evaluated 
numerically  in  terms  of  only  the  mean  p,*  and  the  standard  deviation  a*  of  the  process 
X(t): 


(4=4 
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In  which  Ofu)  is  the  standard  normal  distribution  function. 
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Model  of  Local  Peaks,  Y. 

We  may  instead  choose  to  ignore  all  points  of  the  time  history  except  its  local  peaks, 
typically  defined  as  the  largest  peak  per  upcrossing  of  the  mean  level.  For  a  narrow-band 
normal  process,  this  results  in  a  Rayleigh  distribution  for  Y,  which  again  depends  on  the 
mean  Px  and  the  standard  deviation  Ox  of  the  process  X(t): 


Fy(y)  =  l-exp 


2al 


For  y  >  0  only. 
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Model  of  Global  Peaks,  Z. 

Finally  we  may  instead  choose  the  maximum  value  Z  over  a  still  coarser  time  scale, 
comprising  multiple  peaks  (e.g.,  10-minute  maxima,  1-hour  maxima).  As  when 
proceeding  from  the  process  to  the  local  peaks,  this  step  has  the  advantage  of  focussing 
more  on  the  upper  tail  of  interest,  and  the  corresponding  disadvantage  of  using  less 
detailed  information  about  the  time  history. 

Generally  the  distribution  function  of  Z  is  commonly  estimated  from  that  of  Y  as  follows: 

^z(z)  =  [^j'(y)r  1-4 

In  which  N  here  is  the  number  of  local  peaks  (Y  values)  within  the  duration  over  which  Z 
extends  (again  10  minutes,  1  hour,  etc.),  eq  1.4  assumes  both  that  the  number  of  peaks,  N, 
is  deterministic  and  that  their  levels  are  mutually  independent.  Neither  assumption  is 
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strictly  correct,  but  corrections  commonly  become  less  significant  as  we  consider 
extremes  in  the  upper  tails  of  the  response  probability  distribution.  In  the  Gaussian  case, 
combining  eqs.  1.3  and  1.4  yields  the  results 
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This  finally  gives  us  the  distribution  of  the  extreme  response  based  on  the  mean  and 
standard  deviation  of  the  process  itself  and  the  number  of  local  peaks.  Although  this 
result  is  widely  used,  the  response  of  a  floating  offshore  structure  is  generally  not 
Gaussian  due  to  non-linearities  in  the  response  that  cause  non-Gaussian  behavior  even  if 
the  wave  process  is  Gaussian,  which  it  is  only  by  approximation.  The  various  models  that 
are  available  for  X,  Y  and  Z  in  the  non-Gaussian  case  are  discussed  in  chapter  2.  Given 
arbitrary  distributions  Fyfy)  or  Fz(z)  the  global  extreme  distribution  can  also  be  estimated 
from  the  distribution  of  shorter-period  global  extremes  similar  to  eq.  1.4. 
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If  Fy  has  been  fit  we  use  the  top  expression  in  equation  1.6,  in  which  Ny  is  the  number  of 
local  peaks  expected  in  time  T.  If  Fz  has  been  fit,  we  use  the  bottom  expression  in 


equation  1.6,  in  which  Nz  is  the  number  of  global  peaks  (e.g.  number  of  10-minute  or  1- 
hour  segments)  in  time  T. 


Uncertainty  Estimates  through  Bootstrapping 

Finally,  bootstrapping  methods  (e.g.,  Efron  and  Tibshirani,  1993)  are  used  here  to 
estimate  the  statistical  uncertainty  associated  with  any/all  of  our  estimated  statistics  of 
Xmax-  The  method  is  conceptually  straightforward,  generating  multiple  “equally  likely” 
data  sets  by  simulating,  with  replacement,  from  the  original  data  set.  Thus  some  of  the 
data  values  will  be  repeated  multiple  times,  while  others  will  be  omitted,  in  any  single 
bootstrap  sample  (which  is  of  the  same  size  as  the  original  data  set).  The  same  estimation 
procedure  performed  for  the  original  data  set  is  repeated  for  each  of  the  bootstrapped 
samples,  and  the  net  statistics  on  the  results  are  collected  and  reported. 

The  bootstrap  method  is  “non-parametric”  by  definition,  in  that  it  operates  with  no 
additional  information  beside  the  actual  data  values.  Alternative  approaches  might  fit  a 
parametric  model,  either  statistical  or 

physical,  to  generate  additional  “equally  likely”  samples  from  which  to  infer  sampling 
variability  levels.  Such  approaches  may  confer  advantages  in  some  cases  but  are 
generally  problem-specific;  the  prime  virtue  of  the  bootstrap  method  lies  in  its  generality. 
We  will  consider  both  bootstrapping  and  parametric,  simulation-based  approaches  here. 
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Problem  Description 


Model  Tests 

For  floating  offshore  structures  model  tests  are  often  performed  in  order  to  verify  and 
calibrate  a  computational  model  that  can  reproduce  the  observed  times  trace  of  the 
response  from  the  observed  wave  history.  In  addition  the  results  can  be  used  to  directly 
estimate  the  extremal  statistics,  and  the  distribution  of  the  global  peak  of  the  response  for 
a  certain  sea  state  (i.e.  for  random  realizations  from  a  specific  wave  spectrum 
characterized  by  Hs  and  Tp).  An  important  question  that  arises  for  model  tests  is  how 
accurate  the  results  are,  or  more  precisely  how  much  data  must  be  simulated  to  obtain  the 
required  accuracy. 


Statistics  of  the  Extreme  Response 

As  mentioned  earlier,  the  highest  response  Xmax  in  a  certain  period  of  time  will  be  a 
random  variable  with  a  certain  distribution.  We  may  be  concerned  with  estimating 
various  statistics  of  this  distribution;  e.g.  the  mean  hourly  max  pfXmax.)  its  standard 
deviation  a(Xmax.)  upper  fractiles  of  Xmax.  etc.  Once  such  an  estimate  is  made,  a  measure 
of  the  variability  of  this  estimate,  e.g.  its  standard  deviation,  would  be  a  second  important 
number  to  report.  If ,  for  example,  we  have  10  hours  of  data  for  a  certain  phenomenon 
and  are  interested  in  estimating  the  mean  hourly  maximum,  |J.(Xmax.)  ,  we  would  have  10 
“observed”  hourly  extremes  X1...X10,  from  which  we  can  determine  the  sample  mean, 

X  ,  and  the  sample  standard  deviation,  -xf  .  The  sample 
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mean  would  be  our  best  estimate  of  the  1-hour  mean  max,  li(Xmax»)-  The  standard 
deviation  of  the  mean,  which  is  equal  to  the  standard  deviation,  sx,  divided  by  the  square 
root  of  the  number  of  observations,  would  measure  the  accuracy  of  our  estimate. 

The  first  observation  is  that  we  use  very  little  of  the  available  data,  i.e.  only  one  value  per 
hour.  By  using  very  few  numbers  to  estimate  the  mean  or  any  other  statistic,  we  have  a 
relatively  large  amount  of  statistical  uncertainty.  The  question  now  arises  whether  it  is 
possible  to  use  more  of  the  available  data  by  fitting  a  model  of  e.g.  the  local  peaks,  and 
deriving  the  distribution  of  the  relevant  global  peaks,  in  order  to  reduce  statistical 
uncertainty  of  the  hourly  extreme  in  this  case.  The  price  we  pay  is  that  the  more  data  we 
use,  the  less  relevant  the  included  data  points  will  be  to  the  extreme  value  we  are  trying  to 
estimate,  increasing  model  uncertainty. 

The  Hybrid  Method 

The  hybrid  model  addresses  the  problem  of  data  usage  by  fitting  a  parametric  distribution 
of  the  extreme  response  to  each  hour  in  the  previous  example.  As  discussed  earlier,  this 
approach  allows  any  amount  of  data  between  the  number  of  parameters  and  the  entire 
data  set  to  be  included  in  the  estimate,  by  modeling  the  process,  the  local  peaks,  or  the 
global  peaks.  From  the  fitted  distribution  of  Y  or  Z,  the  distribution  of  the  extreme 
response  and  its  mean  for  that  specific  hour  can  be  determined.  The  hope  is  now  that  the 
variability  in  these  10  estimated  means  is  less  than  the  variability  of  the  10  observed 
extreme  values,  while  they  are  estimates  of  the  same  value.  The  expectation  is  that  there 
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is  an  optimal  combination  between  the  amount  of  data  used  to  fit  the  model  and  the 
relevance  of  the  data  to  the  prediction. 

Torhaug  (1995)  applied  the  hybrid  method  to  determine  the  extreme  response  of  jack-up 
structures,  specifically  modeling  either  all  peaks,  Y,  or  the  entire  process,  X.  In  this  case, 
he  found  that  the  variability  could  be  decreased  by  a  factor  of  2.5,  which  would  reduce 
the  amount  of  required  data  to  achieve  the  same  accuracy  with  a  factor  of  6. 

Scope  and  Organization 

Chapter  2  discusses  methods  to  fit  the  data,  and  probability  distributions  for  the  process, 
the  local  peaks,  and  the  global  peaks.  The  models  that  are  discussed  were  implemented  in 
a  computer  routine  “MaxFits”  (De  Jong  and  Winterstein,  1998).  The  code  is  based  on  a 
distribution  library  that  is  eilso  used  for  the  “Frrs”  routine  (Kashef  and  Winterstein, 

1998).  The  program  was  used  to  analyze  OTRC  spar  buoy  data  to  address  the  problem 
described  in  this  chapter,  and  to  explain  and  demonstrate  the  routine  for  general  use. 

Chapter  3  describes  the  analysis  of  the  spar  buoy  data  in  detail.  It  discusses  the 
distributions  of  the  process,  the  local  peaks,  and  the  global  peaks  respectively.  Estimates 
of  the  uncertainty  are  made,  and  the  various  methods  are  compared  using  1  hour  of 
measured  data.  Chapter  4  compares  the  results  of  chapter  3  to  results  based  on  a  much 
larger  amount  of  simulated  data  in  order  to  reduce  statistical  uncertainty.  The  chapter 
concludes  with  recommendations  for  “optimal”  data  use.  The  final  chapter  summarizes 
the  findings  of  this  work,  and  makes  reconunendations  for  future  research.  It  is  important 


to  note  that  this  work  is  limited  to  the  statistical  analysis  of  stationary  relatively  short 
periods  or  sea  states.  The  results  are  intended  to  be  used  in  long-term  analyses,  for  which 
various  routines  have  been  developed  such  as  “Htcntr”,  “Iform”,  and  “Surfit”,  which 
is  discussed  by  Engebretsen,  1998. 
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Solution  Methods  and  Approaches 


Fitting  distribution  parameters 

All  probability  distributions  used  here  require  a  limited  number  of  parameters  to  be 
determined  in  way  that  allows  us  to  describe  the  data  ‘best’.  What  is  best  is  a  somewhat 
subjective  issue,  although  there  are  several  “goodness  of  fit”  measures  available  (e.g. 
Benjamin  and  Cornell,  1971).  The  following  fitting  methods  can  be  used: 

•  Least  squares  method 

•  Maximum  likelihood  estimate 

•  Method  of  moments 

The  least  squares  fit,  which  returns  distribution  parameters  that  minimize  the  sum  of  the 
squares  of  the  differences  between  the  model  and  the  data,  is  commonly  used.  One  can 
imagine  that  a  small  change  in  one  of  the  parameters  will  hardly  affect  predictions  in  e.g. 
the  body  of  the  distribution,  but  might  have  a  very  significant  effect  on  the  tail  of  the 
distribution.  Turning  this  logic  around,  different  data  points  are  not  equally  “weighted”  in 
the  fit.  Their  “weight”  depends  on  the  type  of  distribution  and  is  completely  arbitrary. 
This  method  is  not  preferred  for  this  reason.  An  alternative  to  the  least  squares  fit  is  the 
maximum  likelihood  estimate.  This  method  maximizes  a  likelihood  function,  which 
expresses  the  likelihood  of  observing  the  data  set  given  the  set  of  parameters  sought.  An 
analytical  solution  to  this  maximization  problem  is  not  always  available,  in  which  case 
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numerical  methods  must  be  used.  The  method  used  here  is  the  method  of  moments.  The 
parameters  are  related  to  the  statistical  moments,  which  represent  all  data  in  a  consistent 
manner.  The  analytical  expressions  for  these  relations  will  be  worked  out  further  on. 

The  moments  are  defined  as: 


=  f{x)dx 

for  i  =  1,2,..  in  the  continuous  case 

2.1a 

for  i  =  1,2,..  in  the  discrete  case 

2.1b 

Where  //,  is  the  moment  and/fxj  is  the  probability  density  function.  Usually  only  the 
first  four  central  moments  are  used,  from  which  the  true  mean  (|i),  standard  deviation  (o), 
skewness  (tta),  and  kurtosis  ((X4)  can  be  estimated  as  follows: 


2.2 


« tr 


2.3 
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2.4 


1  n 

a3=-S 

n  tr 


V 


x.-m^ 


n  r 


n  1=1 


x,-m^ 


V 


V  ^  y 


2.5 


where: 


niv 


as 


:  sample  mean  of  x 
:  sample  skewness  of  x 
:  sample  size 


Sx  :  sample  standard  deviation  of  x 

34  :  sample  kurtosis  of  x 

(sometimesS  is  subtracted,  as  this  is 
04  of  a  normal  distribution) 


These  expressions  are  the  most  conunonly  used,  but  some  may  be  biased  in  certain  cases. 
For  independently  sampled  data,  the  following  corrections  can  be  shown  to  produce 
unbiased  estimates  of  the  quantities  Ox^,  asOx^,  and  cuOx"*  (the  response  cumulants) 
(Fisher,  1928). 


= 


n-1^  n 


= 


(n-lXn-2)-J^  n 


^x(a4-3)  = 


n‘ 


(n-lXn-2Xn-3) 


«  \P  n 


t2 


2.6 


2.7 


2.8 
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Eqs.  2.1-2.%  may  still  produce  bias  in  their  estimates  of  the  unitless  quantities  and  04. 
No  general  results  (independent  of  distribution  type)  are  available  for  tts  and  04  directly. 


Distributions 

The  choice  of  what  distribution  model  to  use  in  what  situation  is  very  much  a  property  of 
the  data.  For  Gaussian  data  the  models  discussed  in  chapter  1  suffice.  In  the  non- 
Gaussian  case  we  need  to  employ  different,  more  general  models.  Various  alternative 
distribution  models  for  the  process,  the  local  peaks,  and  the  global  peaks  are  discussed  in 
the  following  paragraphs. 


The  Process,  X 

A  Gaussian  process  will  plot  as  a  straight  line  on  normal  probability  scale.  Non-Gaussian 
data  will  tend  to  show  curvature  on  this  scale.  One  way  to  add  curvature  to  the  model,  is 
to  distort  a  standard  normal  distribution  such  that  it  plots  as  a  polynomial  on  normal 
probability  scale.  It  has  been  shown  by  Winterstein,  1988)  that  a  cubic  transformation  of 
the  standard  normal  distribution  can  often  capture  non-Gaussian  aspects  due  to  various 
non-linear  mechanisms.  In  this  form  it  is  generally  referred  to  as  the  Hermite  distribution, 
which  can  match  the  first  four  moments  of  the  data.  Another  popular  transformation  of  a 
normal  process  is  to  exponentiate  it,  leading  to  a  “lognormal”  process  (i.e.  one  whose 
logarithm  is  normal).  This  process  can  match  only  two  moments,  however,  and  its  upper 
tail  generally  decays  too  slowly  for  accurate  extreme  estimation. 
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Hermite 


The  Hermite  distribution  is  a  cubic  transformation  of  a  standard  normal  distribution: 

/t/(M)  =  -^exp(-|M^)  2.' 

V2;r 


where: 


x  =  g{u) 

=  m^+KG^\i-Vc^{u^  -l)+C4(w^ ~3m)] 


for  04  >  3 


=  x*-C3(x*-i)-C4(x*'-3jc* 


* 

X 


x-m. 


for  04  <  3 


2.10a 


2.10b 


Approximate  analytical  expressions  are  available  to  relate  C3,  C4  and  k  to  the  skewness 
and  kurtosis  of  the  data  (Winterstein  and  Lange,  1995).  Alternatively  an  iteration  scheme 
can  be  used  to  find  k,  C3,  and  C4. 


Local  Peaks,  Y 

In  the  case  of  a  standard  narrow-band  normal  process  u(t),  a  randomly  selected  peak  v(t) 
will  have  a  Rayleigh  distribution: 


Fv(v)=l-exp(-^v^) 


2.11 


The  Heraiite  model  then  assumes  the  functional  transformation  g  or  g*,  which  relates  u  to 
X  for  each  time  point.  It  therefore  leads  to  the  same  transform  at  the  time  of  peak  values: 


=  ,?(v)  «4  >  3 

v  =  g*(>')  ^4^3 


2.10 


Here  y  describes  a  peak  of  the  non-Gaussian  process,  g  and  g*  are  givem  in  Eq.  10,  and  v 
is  the  Rayleigh  distributed  peak  described  by  Eq.  2.1 1. The  Hermite  model  is  now  used  to 
describe  the  local  peaks,  although  it  has  been  fitted  to  the  process. 

If  we  desire  to  model  the  peaks  directly  and  circumvent  the  distribution  of  the  process 
(using  less  data),  the  Weibull  distribution  is  the  most  commonly  used.  Originally,  the 
Weibull  distribution  was  derived  as  a  type  m  distribution  of  global  minimum  values.  It  is 
based  on  the  assumption  of  peaks  occurring  with  a  limited  tail  of  interest.  Instead  it  is 
often  used  as  a  distribution  of  individual  peaks.  The  Rayleigh  distribution  and  the 
exponential  distribution  are  special  cases  of  the  Weibull. 


Weibull 

The  Weibull  distribution  is  commonly  used  to  model  the  local  peaks  of  a  process. 


Fy(y)=l-exp 


'y-yo 


L  V 


j  j 


2.12 
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If  a  and  P  are  known,  the  mean  and  variance  are  given  by: 


2.14 


These  expressions  can  be  derived  by  transforming  a  standard  exponential  variable  to  a 
Weibull  variable,  and  modifying  the  standard  expressions  for  the  moments  of  a  standard 
exponential  distribution  accordingly.  If  instead  we  estimate  p,y  and  Oy  directly  from  data, 
values  of  a  and  P  have  to  be  determined  from  Eqs.  2.13-2.14  by  iteration. 


The  Weibull  distribution  in  its  form  above  will  plot  as  a  straight  line  on  Weibull-scale.  If 
the  data  show  curvature,  the  Weibull  distribution  can  be  distorted  in  a  way  similar  to  that 
done  for  the  normal  distribution,  leading  to  quadratic  and  cubic  Weibull  distributions. 
The  quadratic  model  allows  the  skewness  of  the  peaks  to  be  fitted  in  addition  to  the  metm 
and  the  standard  deviation.  The  cubic  model  also  allows  the  kurtosis  of  the  peaks  to  be 
fitted.  As  the  distribution  of  peaks  is  one-sided  the  quadratic  Weibull  may  often  suffice. 
By  trying  to  fit  both  the  skewness  and  the  kurtosis  the  result  can  become  unstable  when 
extrapolating  in  the  upper  tail.  This  “overfitting”  of  the  distribution  of  peaks  is  discussed 
in  more  detail  by  Kashef  and  Winterstein  (1998). 

Often  not  all  the  peaks  are  used  to  fit  the  Weibull  distribution.  Instead,  only  the  peaks 
above  a  certain  threshold  are  fitted.  The  threshold,  which  may  often  for  example  reflect 
the  mean  of  the  process,  is  represented  by  the  term  yo  in  expression  2.12.  Alternatively 
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instead  of  adding  a  quadratic  expansion,  this  third  parameter  can  be  fitted  to  the 
skewness,  also  causing  curvature  on  Weibull-scale.  This  however  is  somewhat  indirect  as 
the  lower  threshold  of  the  distribution  is  manipulated  to  better  match  the  upper  tail  of  the 
distribution  to  the  data.  The  quadratic  expansion  on  the  other  hand  manipulates  the  upper 
tail,  and  is  therefore  preferred  here. 


Quadratic  Weibull: 

In  this  model,  one  relates  peaks  y  of  a  process,  possibly  above  some  user-defined  lower 
limit,  yo,  to  a  standard  Weibull  variable  w  by: 


y  -  Jo  =  )  for  tta.y  <  as,, 


for  a3,y  >  a3,v 


2.15 

2.16 


Here  w  has  the  standard  Weibull  distribution: 


Fy(w)  =  l-exp 


f  vv  ' 


2.12 


The  skewness  of  the  model  is  matched  to  the  data  by  changing  e  in  an  iterative  process, 
which  does  not  affect  a  and  p.  These  follow  from  the  previous  expressions  (2.13  and 
2.14). 
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Global  Peaks 


There  are  three  general  types  of  distributions  for  global  extremes;  Type  I,  type  11,  and 
type  ni  distributions.  Their  differences  lie  in  the  assumptions  made  concerning  the 
underlying  distribution  of  peaks.  For  type  I  distributions  we  assume  that  the  upper  tail  of 
the  underlying  distribution  of  Y  falls  off  exponentially,  which  is  e.g.  the  case  for  the 
Weibull  and  the  Rayleigh  distribution.  The  type  I  distribution  is  commonly  called  the 
Gumbel  model,  which  is  described  below. 


Gumbel 

The  Gumbel  distribution  models  the  global  extremes,  Z,  and  has  the  following 
expression: 


(z)  =  exp[-e  ^  ] 


2.18 


where  y  and  Zmode  are  determined  from  the  moments  of  the  global  extremes: 


Tt 


1.282 


2.19  and 


•“modf 


m. 


0.577 

a 


2.20 


where  y  is  Euler’s  constant 

A  Gumbel  distribution  plots  as  a  straight  line  on  Gumbel  probability  scale.  If  the  data 
shows  curvature  a  polynomial  expansion  can  be  added  similar  to  the  Hermite  of  the 
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quadratic  Weibull  distributions.  This  is  discussed  in  detail  by  Winterstein  and  Haver 
(1991). 


Type  II  and  Type  III  distributions 

A  type  n  distribution  assumes  a  limited  tail  on  the  left  at  zero,  but  an  unlimited  tail  on  the 
right  of  the  underlying  distribution  of  individual  peaks.  The  general  form  of  the 
distribution  of  individual  peaks  is: 


Fj.(y)  =  l-y9(y)  fory>0 


2.21 


The  distribution  of  the  global  peaks,  Z,  is  the  distribution  of  Y  raised  to  the  power  of  the 
number  of  individual  peaks,  assuming  independence.  This  gives  after  some 
approximations: 


for  z  >  0 


2.22 


Type  II  distributions  are  related  to  type  I  distributions  in  a  similar  fashion  as  the 
lognormal  distribution  is  related  to  the  normal  distribution,  i.e.  Zj  =  InZn,  where  Zi  has  a 
type  I,  and  Zu  has  a  type  II  distribution.  Contrary  to  Type  I  distributions.  Type  11  and  HI 
distributions  are  ill-suited  to  our  moment-fits,  due  to  potential  moment  divergence  (Type 
11)  or  to  the  difficulty  in  predicting  truncated  distributions  (Type  III)  from  moment 
information. 
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Spar  Buoy  Results  Based  on  Ihr  of  Model  Test  Data 


Introduction 

Recently  the  first  two  production  spars  have  been  built.  The  Neptune,  operated  by  Orja 
was  installed  in  1996,  while  Chevron  installed  its  Genesis  platform  in  1997.  A  third  spar 
is  under  construction  for  the  Diana  field,  which  will  be  operated  by  Exxon.  Presently 
production  spars  are  considered  by  many  operators  as  a  very  competitive  concept  for 
water  depths  ranging  from  roughly  2000  to  6000  ft,  offering  important  advantages  over 
tension  leg  platforms  (TLPs)  and  semi-submersibles. 

In  principle  a  spar  buoy  is  a  very  large  vertical  cylinder,  providing  buoyancy  to  support 
production  facilities,  and  optionally  to  store  oil.  The  structure  is  kept  in  place  by  mooring 
lines,  which  can  be  either  catenary  or  taut.  The  latter  concept  may  require  the  use  of 
polyester  lines  and  is  an  upcoming  solution  for  deep  water.  It  provides  restoring  forces 
through  the  elastic  properties  of  the  material,  in  addition  to  the  restoring  forces  resulting 
from  the  geometry  and  weight  of  the  lines,  which  give  a  catenary  mooring  system  its 
stiffness.  The  horizontal  stiffness  of  the  spar  system  is  determined  by  the  mooring 
system,  while  the  vertical  stiffness  is  determined  by  hydrostatic  properties,  which  are 
determined  by  the  geometry  of  the  hull.  The  moorings  and  the  structure  are  designed  such 
that  natural  periods  are  well  above  the  dominant  wave  periods.  Natural  periods  for  the 
different  modes  of  a  production  spar  are  generally  in  the  order  of  minutes,  while  wave 
periods  are  in  the  order  of  seconds.  Nevertheless  dynamic  excitation  can  occur,  which  is 
discussed  in  the  following  paragraphs.  Damping  is  caused  by  loss  of  energy  through 


turbulence  and  friction,  and  is  generally  difficult  to  quantify.  We  include  here  the  effects 
of  viscous  damping,  radiation  damping  and  wave  drift  damping.  The  mass  of  the  system 
is  the  mass  of  the  spar  buoy  plus  the  mass  of  the  water  entrapped  in  the  center  well,  and 
the  added  mass  of  the  water  around  the  structure.  The  reader  can  refer  to  e.g.  Glanville  et 
al  (1991),  for  more  general  information  about  the  spar  buoy  concept. 


The  Model  Test 

The  data  that  is  analyzed  here  is  the  result  of  a  joint  industry  project,  carried  out  by  the 
Offshore  Technology  Research  Center  at  Texas  A&M.  Model  tests  were  run  for  a  spar 
buoy  with  a  catenary  mooring  system.  The  characteristics  of  the  structure  are  sununarized 
in  table  3.1.  The  prototype  scale,  on  which  all  values  are  reported,  is  1:55. 

The  results  have  been  used  to  study  the  accuracy  of  various  computational  models  that 
can  predict  structural  response  accurately  (e.g.  Ran  et  al.,  1996,  Weggel  and  Roesset, 
1996,  or  Mekha  and  Roesset,  1995).  Issues  that  were  of  particular  interest  were  the 
damping  and  viscous  forces.  Different  conditions  were  simulated,  of  which  one 
realization  of  alOO  year  extreme  unidirectional  Gulf  of  Mexico  (GOM)  storm  without 
current  was  selected  for  the  analysis  here.  The  characteristics  of  this  sea  state  are  listed  in 
table  3.2.  The  following  data  has  been  used  from  these  tests: 

•  The  undisturbed  (far  field)  wave  elevation,  measured  at  an  along-wave  point,  during 
the  actual  model  test. 

•  The  disturbed  wave  elevation  near  the  structure 

•  Horizontal  offset  of  the  spar  at  54.8  m  above  MWL,  along  the  vertical  cylinder  axis 
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•  Pitch  rotation 


•  Heave  motion 

The  horizontal  offset  is  the  combined  effect  from  surge  and  pitch  motions.  In  general 
besides  the  motions  mentioned  above,  also  lateral  horizontal  motion  or  sway,  caused  by 
vortex  induced  vibrations  (VIV)  can  be  significant.  In  practice,  these  motions  are  reduced 
by  helical  strakes  on  the  hull,  which  can  be  seen  in  figure  3.2.  Sway  will  not  be  discussed 
here.  For  the  analysis  only  the  horizontal  offset  data  at  54.8  m  MWL  was  used. 


Description 

Notation 

Value 

Diameter 

D 

40.5  m 

Draft  from  MWL 

H 

198.2  m 

Mass  with  entrapped  water 

M 

2.59- 10^  kg 

MWL  to  center  of  gravity 

Zcg 

105.8  m 

MWL  to  center  of  buoyancy 

Zcb 

99.1  m 

MWL  to  mooring  connection 

Zf 

105.6  m 

Radius  of  gyration  to  MWL 

Kr 

122.8  m 

Horizontal  mooring  stiffness 

K 

191  kN/m 

Location  of  gauge  above  MWL 

Zm 

54.8  m 

Water  depth 

D 

922  m 

Table  3.1  Characteristics  of  test  spar  buoy 
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105.6  m 


Fig  3.1  Cross  section  of  a  spar  buoy 


Fig  3.2  Installation  of  the  Neptune  spar  buoy 


Test 

GOMl 

Npts 

9702 

Nominal  Hs  (m) 

13.1 

Nominal  Tp  (s) 

14.0 

Observed  Hs  (m) 

14.2 

Observed  Tp  (s) 

14.1 

Calculated  Tz  (s) 

10.8 

Mean  (m) 

0.06 

Standard  deviation  (m) 

3.54 

Skewness 

0.31 

Kurtosis 

3.06 

Minimum  (m) 

-9.05 

Maximum  (m) 

12.74 

Table  3.2  Characteristics  of  test  sea  state 
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The  Data 


Fourier  Spectrum 

The  Fourier  spectrum  of  the  total  response  (figure  3.3)  shows  3  distinct  peaks,  at  0.003, 
0.015,  and  0.07  Hz.  These  peaks  correspond  respectively  to  the  natural  periods  in  surge 
and  pitch,  and  to  the  spectral  peak  period  of  the  sea  state  (330s,  67s,  14s).  All  the  wave 
energy  is  at  the  frequencies  around  the  spectral  peak  period,  where  the  response  is  related 
approximately  linearly  to  the  waves.  Although  there  is  no  wave  energy  at  the  two  low 
frequency  peaks  for  surge  and  pitch,  the  system  is  excited  at  its  natural  frequencies  by 
difference  frequencies  in  the  waves.  This  is  the  result  of  a  non-linear  forcing  mechanism, 
which  is  discussed  in  detail  by  e.g.  Faltinsen  et  al.  (1990). 
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Time  Series 


In  order  to  analyze  the  data,  the  response  time  series  was  filtered  into  3  different 
components,  surge,  pitch,  and  wave  frequency  response,  by  simply  introducing  frequency 
cutoffs  at  0.006  and  0.03  Hz.  The  filtered  time  series  for  the  3  components  and  the  total 
are  plotted  in  figure  3.4.  For  the  selected  location  of  the  observed  motions,  the  3 
components  are  of  almost  equal  magnitude  (7-9  m).  Also  note  that  in  surge  there  are  on 
average  only  1 1  peaks  per  hour,  while  in  pitch  there  will  be  54  peaks  on  average.  The 
limited  amount  of  data  of  the  extremes  will  affect  the  predictions  we  intend  to  make 
further  on  in  this  report.  There  have  been  discussions  suggesting  transient  effects  or 
“mode  swapping”  in  the  time  trace  of  the  total  response.  The  observed  behavior  in  this 
test,  however,  was  found  to  be  statistically  consistent  with  one-hour  segments  of 
simulated  steady  state  response  to  stationary  random  seas. 

Moments 

The  moments  of  the  process  and  the  peaks  are  given  in  table  3.3.  The  mean  of  the 
components  by  definition  is  zero,  due  to  the  filtering  operation.  The  mean  total  offset  is 
almost  5  m,  and  has  been  removed  for  all  the  analyses.  Also  note  that  in  particular  for 
surge  peaks,  these  moments  are  estimated  from  very  few  numbers,  and  are  therefore 
“noisy”.  The  limited  amount  of  data  can  also  bias  the  higher  moments  downward. 

The  base  case  we  will  refer  to  is  a  Gaussian  process,  for  which  the  probability 
distribution  of  the  continuous  time  trace  will  be  normal.  For  a  normal  distribution  the 

f 

mean  and  standard  deviation  determine  the  distribution. 
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The  higher  moments,  the  skewness  and  kurtosis,  are  implicit  in  the  normal  distribution. 
Their  values  are  respectively  0  and  3.  In  contrast,  the  columns  denoted  m3  and  m4  in 
table  3.3  report  sample  skewness  and  kurtosis  values,  estimated  as  in  Eqs.  2.7-2.8.  From 
these  values,  one  may  be  tempted  to  infer  that  the  responses  show  narrower  than 
Gaussian  distributions:  although  the  skewness  values  of  the  process,  reported  in  table  3.3, 
are  nearly  0,  the  kurtosis  values  are  systematically  less  than  3.  We  should  note,  however, 
that  these  low  values  of  the  kurtosis  may  be  due,  at  least  in  part,  to  the  effect  of  limited 
data.  This  lack  of  data  will  commonly  lead  to  systematic  underestimation  of  kurtosis 
values.  Note,  also,  that  the  smallest  estimated  kurtosis  values  (2.17  and  2.20)  are 
associated  with  pitch  and  surge,  for  which  fewer  independent  cycles  of  motion  are 
available. 


Process 

mean 

stdv 

m3 

m4 

Surge 

0.00 

3,37 

0.10 

2.20 

Pitch 

0.00 

4.00 

-0.01 

2.17 

Wavefr. 

0.00 

2.47 

-0.03 

2.70 

Total 

4.94 

5.79 

0.08 

2.68 

Peaks 

mean 

stdv 

m3 

CoV 

Surge 

4.39 

1.93 

-0.14 

0.44 

Pitch 

5.26 

2.02 

-0.03 

0.38 

Wavefr. 

3.11 

1.55 

0.29 

0.50 

Total 

5.71 

4.41 

0.76 

0.77 

Table  3.3  Moments  of  the  response  processes  and  their  local  peaks 


Alternatively,  of  course,  these  apparent  narrower-than-Gaussian  statistics  may  be  due  to 
some  physical  hardening  mechanism,  such  as  wave  drift  damping  or  perhaps  non¬ 
linearity  in  the  mooring  forces.  These  two  hypotheses  (statistically  biased  moments  due 
to  limited  data,  vs.  true  hardening  behavior)  will  be  considered  further  in  chapter  4, 
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which  uses  long  simulated  histories  from  a  theoretical  spar  buoy  response  model  to 
quantify  the  bias  due  to  limited  data. 


As  noted  in  chapter  2,  if  the  time  trace  of  the  process  is  a  narrow  banded  Gaussian 
process,  its  peaks  will  follow  a  Rayleigh  distribution.  The  Rayleigh  distribution  of  peaks 
corresponsing  to  a  general  normal  process  is  given  below,  where  jUx  represents  the  mean 
of  the  process,  and  Ox  the  standard  deviation  of  the  process. 


The  Rayleigh  distribution  represents  a  special  case  of  the  Weibull  distribution  (Eq.  2.12), 
where  the  exponent,  a  is  equal  to  2,  P  is  equ£il  to  axV2,  and  y  =  Px-  The  coefficient  of 
variation  of  the  peaks  implicit  in  the  Rayleigh  distribution  when  px  =  0  is  equal  to  0.5228 
(Johnson  &  Kotz,  1970).  The  implicit  skewness  is  0.631 1  (Johnson  &  Kotz,  1970).  As 
noted  in  chapter  2,  the  kurtosis  of  the  peaks  contains  little  additional  information  here,  as 
it  is  most  useful  for  two-sided  distributions,  while  the  distribution  of  peaks  is  only  one¬ 
sided  (Kashef,  1998). 

Table  3.3  shows  that  a  Rayleigh  model  of  peaks  is  most  nearly  consistent,  in  its  CoV 
value  (0.52),  with  the  wave  frequency  response  peaks.  As  the  wave  frequency  response  is 
physically  thought  to  be  a  nearly  linear  transformation  of  the  nearly  Gaussian  wave  input, 
this  result  is  not  surprising.  In  contrast,  surge  and  pitch  results  suggest  somewhat 
narrower  distributions,  in  terms  of  slightly  lower  CoV  values  and  notably  smaller 


skewnesses  than  the  Rayleigh  model.  Again,  however,  some  of  these  effects  may  be  due 
to  the  limited  model  test  duration  of  1  hour,  which  most  severely  impacts  the  statistics  of 
surge  and  pitch  components.  In  contrast,  the  peaks  of  the  total  process  show  a 
significantly  higher  CoV  (0.77)  than  any  of  the  components.  As  will  be  discussed  further 
below,  this  arises  due  to  the  wide-band  nature  of  the  total  response:  its  mixture  of  3 
frequency  components  leads,  roughly  speaking,  to  a  mixture  of  Weibull  distributions, 
with  correspondingly  larger  CoVs. 

The  foregoing  results  at  least  suggest  the  need  to  include  the  potential  for  (1)  non- 
Gaussian  distributions  of  the  process  components;  and  (2)  non-Rayleigh  distributions  of 
the  response  peaks.  These  are  investigated  in  the  next  section,  using  both  non-Gaussian 
(Hermite)  process  models  and  non-Rayleigh  (Weibull  and  quadratic-Weibull)  models  of 
peaks. 


Distribution  of  the  Peaks 

The  4  distribution  models  that  were  discussed  in  the  previous  chapter  were  fitted  to  the 
spar  buoy  data  in  order  to  study  their  applicability  before  predicting  extreme  values, 
which  is  our  objective.  The  Hermite  transformed  Rayleigh  and  both  Weibull  distributions 
describe  the  distribution  of  local  peaks  and  can  be  compared  directly,  although  the 
Hermite  modelis  fitted  to  the  moments  of  the  original  process.  The  Gumbel  distribution 
describes  the  global  peaks  and  therefore  cannot  be  compared  directly  to  the  other 
distributions  of  local  peaks.  The  number  of  global  peaks  is  selected  by  the  user,  and 
requires  some  insight  in  the  data.  In  this  case  surge  is  the  governing  factor  with  the  least 
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amount  of  peaks.  Selecting  more  global  peaks  than  there  are  local  peaks  would  be 
inconsistent,  as  not  all  of  the  values  in  the  set  of  ‘selected  peaks’  will  be  real  or  local 
peaks.  We  selected  the  10-minute  global  peaks  for  all  time  series,  which  gives  6  jieaks 
per  hour  as  a  compromise  between  obtaining  meaningful  global  peaks  and  a  sufficient 
amount  of  data.  The  results  were  somewhat  sensitive  to  the  amount  of  global  peaks.  It 
turns  out  however  that  this  effect  disappears  if  we  have  multiple  hours  of  data  (see 
section  4.4). 


Distribution  Plots 

The  slope  of  the  local  peak  data,  when  plotted  on  Weibull,  scale  is  equal  to  the  exponent 
a.  The  coefficient  of  variation  can  be  estimated  roughly  by  1/a.  The  slope  of  data  on 
Weibull-scale  gives  us  an  impression  of  how  nearly  the  data  follow  a  Rayleigh  model 
(a~2).  The  slope  of  the  data  on  Gumbel-scale  is  equal  to  -a,  which  controls  the  standard 
deviation  of  the  Gumbel  distribution  in  an  inverse  way.  In  particular,  the  standard 
deviation  of  the  Gumbel  model  equals  1.282/a.  Curvature  on  either  scale  suggests  the 
effect  of  higher  moments  is  important  to  capture. 

Local  Peaks 

The  distributions  of  local  peaks  are  shown  in  Figures  3.5-3. 8.  A  first  observation  is  that  if 
we  exclude  the  case  of  the  total  response  (Figure  3.8),  local  peak  data  in  the  remaining 
three  cases  show  no  significant  curvature  on  these  Weibull  scale  plots.  Hence,  a  Weibull 
model  appears  to  suffice  for  all  three  components  —surge,  pitch,  and  wave  frequency 
response.  The  higher-moment  models  (quadratic  Weibull,  Hermite)  lead  to  no  significant 
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difference,  and  do  not  appear  in  these  cases  to  be  worth  the  additional  parameter(s)  they 
require  to  fit.  Among  these  three  components,  the  Weibull  models  of  both  surge  and 
wave-frequency  response  have  slopes  a  near  2  (see  Table  3.4),  hence  a  still  simpler 
Rayleigh  model  may  suffice  in  these  cases.  The  pitch  component  shows  a  narrower-than- 
Rayleigh  distribution,  with  alpha  nearly  doubled  (3.7),  so  that  the  best  Weibull  model 
provides  a  notable  improvement  over  the  Rayleigh  in  this  case. 


Weibull 

mean 

stdv 

a 

p 

Surge 

4.27 

2.16 

^07 

4.82 

Pitch 

5.80 

1.75 

3.69 

6.43 

Total 

5.99 

4.67 

1.29 

6.48 

Wawefr, 

3.09 

1.57 

3.49 

Q  Weibull 

MTin 

1C 

1/a 

e 

Surge 

-1.34 

8.88 

1.13 

0.48 

-0.81 

Pitch 

-0.29 

6.27 

1.11 

0.27 

-0.03 

Total 

-1.85 

12.66 

1.08 

0.77 

-0.73 

Wawefr. 

-0.55 

4.59 

1.13 

0.48 

-0.28 

Hennite 

mean 

stdv*K 

c3 

C4 

Surge 

-0.11 

3.26 

0,03 

-0.06 

Pitch 

0.01 

3.99 

-0,01 

-0.12 

Total 

4.84 

5,94 

0.02 

-0.02 

Wawefr. 

0.00 

2.46 

0.00 

-0.02 

Table  3.4  Distribution  parameters 


Finally,  we  consider  the  total  response  peaks  (Figure  3.8),  which  are  the  quantities  of 
ultimate  interest.  Here  the  data  itself  does  not  convincingly  follow  a  single  Weibull 
model  (i.e.,  appear  to  plotiinearly)-throughout  its  entire  range.  It  might  appear  that  a  tail- 
fit  Weibull  model  (e.g.,  to  above  9  m)  would  lead  to  a  somewhat  steeper  line  (narrower 
distribution)  than  the  moment-fit  Weibull  model  shown  in  the  figure.  To  explain  this,  it  is 
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useful  to  recall  that  the  Rayleigh  model,  and  its  Weibull  generalizations,  were  postulated 
originally  as  models  of  peaks  of  relatively  narrow-band  processes. 

We  may  expect  this  narrow-band  assumption  to  be  fairly  accurate  for  each  of  the  three 
component  processes.  The  total  response,  however,  represents  a  mixture  of  these  three 
components,  notable  in  view  of  (1)  their  very  different  time  scales  and  (2)  their 
similar  contributions  to  the  total  variance  of  the  process  (Table  3.3).  The  resulting 
process  will  be  far  from  narrow-banded,  and  its  peaks  may,  roughly  speaking,  be 
considered  as  some  type  of  mixture  of  Weibull  distributions.  In  such  cases,  it  may  be 
difficult  to  expect  a  single,  moment-fit  Weibull  to  accurately  track  the  distribution  of 
peaks  far  into  its  upper  tail.  The  added  flexibility  of  the  quadratic  Weibull  may  prove 
useful  in  this  case.  Note  too  that  the  Hermite  model  can  be  argued  to  be  an  adequate 
model  in  the  upper  part  of  the  distribution,  though  it  far  overestimates  peaks  at  lower 
fractile  levels.  This  is  not  accidental.  Recall  that  the  Hermite  model  does  not  sample  the 
actual  peaks  at  all;  it  merely  samples  moments  of  the  non-Gaussian  process,  and  creates  a 
consistent  narrow-band  process  by  transforming  a  Rayleigh  model  of  peaks.  Because  the 
total  response  is  in  fact  wide-band,  the  narrow-band  Hermite  model  cannot  (and  does  not 
seek  to)  follow  the  peak  distribution  throughout  its  entire  body.  Because  the  effect  of 
bandwidth  becomes  lessened  for  high  response  peaks,  however,  it  may  still  be  hoped  that 
the  Hermite  model  remains  adequate  for  rare  peak  levels  of  practical  interest.  (This  will 
again  be  taken  up  in  Chapter  4,  where  additional  data  better  illustrate  the  differences 
among  these  distributions  of  total  response  peaks.) 
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10-minute  Global  Peaks 


Finally,  Figures  3.9  and  3.10  show  the  observed  distribution  of  global  peaks:  i.e.,  the  6 
10-minute  maxima  in  the  1-hour  history.  These  figures  are  plotted  on  Gumbel  scale;  also 
shown  in  each  case  is  the  straight  line  arising  from  a  Gumbel  fit  to  the  mean  and  standard 
deviation  of  these  6  maxima.  With  6  data  here,  there  is  little  chance  to  assess  deviations 
from  the  Gumbel  model.  What  is  clear,  however,  is  that  the  10-minute  maxima  from  the 
wave  frequency  response  (upper  half  of  Fig  3.10)  show  considerably  less  variability 
(higher  slope)  than  all  other  cases. 

In  other  words,  once  we  compare  global  maximum  responses  over  the  same  interval 
(here,  10  minutes),  quantities  involving  slow  drift  motions  (pitch,  surge,  and  total 
response)  show  greater  variability  than  those  due  to  wave-frequency  motions.  This  is 
due  in  part  to  the  far  fewer  cycles  of  slow  drift  motion,  and  also  to  the  light  damping 
which  induces  high  correlation  among  peaks  (see  table  3.5).  Note,  for  example,  the 
episodic,  slowly  evolving  nature  of  the  surge  component  of  this  model  test,  as  shown  in 
Figure  3.4.  It  is  not  difficult  to  imagine  that  its  10-minute  extremes  will  show  far  more 
scatter  than  those,  for  example,  of  the  wave-frequency  response  shown  in  the  same 
figure.  The  net  impact  of  this,  as  we  shall  see  in  following  sections,  is  that  we  will 
generally  require  longer  data  sets  to  accurately  estimate  extremes  of  slow-drift  motions, 
as  compared  with  cases  of  wave-frequency  motion  alone. 


Peaks: 

global 

local 

Surge 

0.33 

0.50 

Pitch 

-0.07 

0.69 

Total 

-0.13 

0.03 

Wave  fr. 

0.00 

0.35 

Table  3.5  Correlation  coefficients  of  subsequent  local  peaks 


35 


probability  of 


1 


2  3  4  56789  10 

offset  (m) 


Fig  3.5  Surge  data  and  fitted  distributions  plotted  on  WeibulUscale 
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probability  of  exceedence 


Pitch  (1hr) 


offset  (m) 


Fig  3.6  Pitch  data  and  fitted  distributions  plotted  on  WeibuU-scale 
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probability  of  exceedence 


Wave  Frequency  Response  (1hr) 


1  2  3  456789  10 

offset  (m) 


Fig  3.7  V\Iave  frequency  response  data  and  fitted  distributions  plotted  on  WeibulUscale 
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probability  of  exceedence 


Total  Response  (Ihr) 


Fig  3.8  Total  response  data  and  fitted  distributions  plotted  on  Weibull-scale 
(Note  different  horizontal  scale  from  components) 
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probability  of  exceedence  probability  of  exceedence 


Pitch  (Ihr) 


offset  (m) 

Fig  3.9  Data  and  fitted  distributions  plotted  on  GumbeUscale  for  surge  and  pitch 
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probability  of  exceedence  probability  of  exceedence 


Wave  Frequency  Response  (1hr) 


offset  (m) 

Fig  3.10  Data  and  fitted  distributions  plotted  on  Gumbel-scale  for  the  total  and  wave  frequency  response 
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Estimating  Statistics  of  the  Extreme  Value  Distribution 

This  section  discusses  how  we  use  the  foregoing  results  to  estimate  the  quantity  of 
ultimate  interest;  namely,  statistics  of  the  distribution  of  Xmax,  the  extreme  response  to 
occur  in  an  arbitrary  1-hour  duration. 

As  noted  earlier,  a  first  quantity  of  interest  is  lixmax,  the  average  value  of  Xmax  to  be 
expected  in  an  arbitrary  hour.  From  the  single  hour  of  model  tests,  p,xmax  can  be  directly 
estimated  by  the  observed  maximum  value  in  that  hour.  (This  will  of  course  be  quite  a 
"noisy"  estimate,  but  it  will  at  least  be  unbiased.)  Table  3.6  shows  these  direct  estimates 
of  the  mean  max,  for  each  component  and  the  total  response,  in  the  row  marked 
"Observed".  Other  entries  in  the  columns  "Mean  max"  give  estimates  based  on  fitting  one 
of  the  foregoing  models  (Hermite  of  process,  Weibull  or  quadratic  Weibull  of  local 
peaks,  Gumbel  model  of  10-minute  peaks)  to  the  hour  of  data.  To  do  this,  note  that  each 
model  first  estimates  the  distribution  function  of  Xmax,  based  on  either  the  distribution 
function  of  Y  or  Z  (as  in  Eq.  1.6).  The  mean  value  of  Xmax  is  then  estimated  from  its 
distribution  function  by  numerical  quadrature  methods. 
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Surge 


Wave  Frequency  Response 


GOMl 

Mean  max 

Bs^Stdv 

Cov 

#peak5 

GOMl 

Mean  max 

Bs^Stdv 

OoV 

#peak5 

Obsened 

7^4 

0.46 

0.064 

11 

Obsen^ 

7.16 

0.17 

0.024 

265 

GLfTt)el 

7.93 

0.89 

0.113 

6 

Gumbel 

7.31 

0.29 

0.040 

6 

Weibiii 

7.88 

0.77 

0.097 

11 

Weibdl 

a32 

0.30 

0.036 

265 

Q  Weibull 

7.54 

0.65 

0.086 

11 

QWatxil 

7.66 

0.27 

0.035 

265 

Herrnte 

7.16 

- 

- 

11 

Hemrite 

7.62 

* 

• 

265 

Pitch 

1  Total  Response  | 

GOMl 

Mean  max 

Bs_Stdv 

CoV 

#peaks 

GOMl 

MeannrBx 

Bs^Stdv 

CdV 

#peek5 

Observed 

9.52 

0.33 

0.035 

54 

Observed 

17.82 

0.72 

0.040 

112 

Gintei 

10.03 

0.99 

0.059 

6 

Gumbel 

17.86 

0.91 

0.051 

6 

Weibull 

10.12 

0.56 

0.054 

54 

WeifcxJI 

22.30 

1.41 

0.063 

112 

Q  Weibull 

9.59 

0.46 

0.048 

54 

QWatxil 

19.49 

1.10 

0.057 

112 

Hemrite 

9.31 

• 

• 

54 

Hernite 

17.34 

- 

- 

112 

Table  3.6  Mean  hourly  extremes  for  GOMl 


Bias 

In  principle,  the  bias  associated  with  any  model’s  estimate  of  |Xxmax  can  be  estimated  from 
the  ratio  between  its  numerical  estimate  and  the  "Observed"  unbiased  value.  As  our 
observed  value  uses  only  a  single  observed  extreme,  however,  we  should  not  interpret 
these  bias  estimates  too  literally.  It  is  somewhat  encouraging,  though,  that  most  estimates 
lie  within  10%  of  the  observed  values.  The  notable  exception  is  the  Weibull  estimate  for 
the  total  response,  which  exceeds  the  observed  value  by  25%.  This  perhaps  supports  the 
previous  observation  that  since  the  total  response  shows  three  distinct  spectral  modes,  its 
peaks  are  less  well  modeled  by  a  single  moment-fit  Weibull  model  (as  compared  with  the 
peaks  of  the  relatively  narrow-band  component  processes).  This  will  be  borne  out  in 
Chapter  4,  in  which  the  only  significant  bias  results  from  applying  the  Weibull  model  to 
estimate  extremes  of  the  total  response. 
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Uncertainty 

In  addition  to  seeking  an  estimate  that  is  unbiased  —  i.e.,  correct  on  average  --  we  also 
wish  to  consider  the  variability  in  the  estimate  due  to  limited  data.  In  particular,  we  wish 
to  determine  whether  fitting  a  model  is  useful  in  reducing  variability  in  the  estimated 
mean  maximum  {ixmax.  as  compared  with  using  simply  the  observed  maximum  value. 

As  noted  in  Chapter  1,  bootstrapping  provides  one  method  of  estimating  this  variability. 
Its  central  idea  is  to  create  many  synthetic,  "equally  likely"  one-hour  samples,  by 
sampling  from  the  observed  data  set  with  replacement.  For  example,  the  1-hour  surge 
component  response  is  found  here  to  have  1 1  local  peaks  (i.e.,  peaks  above  the  mean 
response  level).  As  shown  in  Table  3.6,  a  Weibull  model  applied  to  these  1 1  data  yields 
the  estimate  of  |Xxmax  =  7.88m.  Bootstrapping  methods  then  suggest  that  this  estimate  has 
a  standard  deviation  of  0.77m,  and  hence  a  CoV  of  .097  (the  columns  marked  "Bs_Stdv" 
and  "CoV"  in  Table  3.6). 

To  understand  the  bootstrap  process  used  here,  imagine  that  we  have  an  1 1 -sided  die, 
with  each  of  the  1 1  local  peak  values  written  on  one  of  its  faces.  We  then  imagine 
independently  rolling  this  die  1 1  times.  The  resulting  1 1  outcomes  are  regarded  as  a  new 
"bootstrap"  sample.  (Each  of  the  outcomes  must  be  a  value  of  the  original  sample; 
moreover,  it  is  almost  inevitable  that  some  values  will  arise  multiple  times,  and  other 
values  omitted,  from  any  particular  sample.)  We  then  repeat  this  exercise  to  create  many 
bootstrap  samples,  and  keep  track  of  the  estimates  of  p-xmax  we  find  from  each  sample. 
The  standard  deviation  of  these  estimates  (here  based  on  200  bootstrap  samples) 
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has  been  reported  in  Table  3.6.  The  quadratic  Weibull  model  of  surge  is  also  based  on  all 
1 1  of  its  peaks,  hence  its  standard  deviation  is  also  found  by  creating  new  bootstrap 
samples  of  these  1 1  values.  In  contrast,  the  Gumbel  model  is  based  on  an  observed 
data  set  of  only  6  10-minute  maxima;  these  form  the  basis  of  new  bootstrap  samples  (of 
size  6),  to  which  the  Gumbel  estimator  is  repeatedly  applied  to  assess  its  variability. 

(The  standard  deviation  associated  with  the  observed  values  of  P-xmax  are  found  here  by 
bootstrapping  the  entire  vector  of  local  peaks,  again  1 1  in  the  case  of  surge,  and  then 
selecting  the  maximum  of  each  bootstrap  sample.) 

In  surveying  the  bootstrap-estimated  CoV  values  in  Table  3.6,  one  may  note  first  that  the 
observed  estimates  apparently  show  less  variability  than  those  produced  by  any  of  the 
corresponding  models.  If  these  results  are  truly  representative,  they  would  suggest  that 
fitting  these  models  is  in  fact  disadvantageous.  But  in  addition,  the  net  CoV  values  seem 
unusually  low,  particularly  those  associated  with  the  direct,  observed  estimates  (e.g.,  only 
2.4%  in  wave-frequency  extreme,  and  4%  in  pitch  and  total  extremes).  One  may  question 
whether  these  numbers  are  truly  indicative  of  long-run  behavior,  or  perhaps  unusually 
low  due  to  (1)  the  happenstance  of  the  particular  hour  history  chosen;  or  (2)  some 
systematic  uncertainty  underestimation  by  the  bootstrap  method  as  currently 
implemented.  Regarding  point  (2),  our  simple  bootstrap  sampling  of  response  peaks 
neglects  their  correlation,  which  may  well  serve  to  underestimate  the  variability  due  to 
limited-duration  samples.  In  addition,  the  bootstrapping  method  is  in  general  tail-limited, 
failing  to  sample  beyond  the  last  data  point  observed,  which  may  also  affect  extreme 
value  estimation. 
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As  in  other  issues  raised  in  this  chapter,  distinction  between  issues  (1)  and  (2)  above 
requires  a  longer  database  of  results,  from  which  systematic  trends  can  be  distinguished 
from  the  particular  values  of  this  limited  1-hour  test.  This  is  the  topic  of  Chapter  4. 
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probability  of  exceedence  probability  of  exceedence 


1 


Distribution  of  the  Hourly  Extreme  Surge  Response 


Distribution  of  the  Hourly  Extreme  Pitch  Response 


offset  (m) 

Fig  3.11  Distribution  of  the  predicted  hourly  extreme  value  for  surge  and  pitch 
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probability  of  exceedence  probability  of  exceedence 


Distribution  of  the  Hourly  Extreme  Wave  Frequency  Response 


Distribution  of  the  Hourly  Extreme  Total  Response 


offset  (m) 

Fig  3.12  Distribution  of  the  predicted  hourly  extreme  value  for  the  total  and  wave  frequency  response 
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Spar  Buoy  Results  Based  on  36  hrs  of  Data 


Introduction 

The  results  in  the  previous  chapter  are  subject  to  a  great  deal  of  statistical  uncertainty  due 
to  the  very  limited  amount  of  suitable  model  test  data  that  was  available.  The  two  main 
issues  that  emerged  from  chapter  3:  “Are  bootstrap  estimates  of  the  extreme  statistics 
accurate?”  and  “Are  the  predictions  of  the  distribution  of  extreme  values  made  in  chapter 
3  biased?”  will  be  revisited  using  an  amount  of  data  that  intends  to  eliminate  this 
uncertainty.  In  addition  we  will  explore  the  answer  to  our  original  question  “Does  fitting 
a  model,  utilizing  more  of  the  original  data,  reduce  the  uncertainty  in  our  predictions  of 
the  extreme  values?”  further. 

In  order  to  evaluate  the  results  bootstrapping  provided  in  the  previous  section,  we  seek  to 
create  multiple  data  sets  by  numerically  simulating  the  model  test,  instead  of  resampling 
from  the  original  data  set.  Although  the  model  test  that  was  analyzed  was  repeated  once, 
the  two  hours  of  measured  data  do  not  provide  a  sufficient  amount  of  data  to  improve  the 
results  of  our  analyses  much.  For  the  computer  simulation  we  used  the  results  from 
second  order  diffraction  analysis,  and  post-processed  these  using  a  routine  called 
“TFPOP”  (Ude  1996).  In  fact  this,simulation-based  method  may  be  regarded  as  a  form  of 
parametric  bootstrapping  as  opposed  to  the  non-parametric  bootstrapping  method  used  in 
the  previous  sections.  The  TFPOP-model  for  the  OTRC  spar  buoy  was  developed  and 
verified  by  Jha  ( 1 997). 
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Non-parametric  bootstrapping  may  not  be  very  useful  to  estimate  the  bias  in  our 
predictions  of  the  mean  extreme  value,  as  all  the  data  sets  we  compare  are  produced  with 
the  same,  possibly  biased  estimator.  Simulation  however,  does  allow  us  to  evaluate  the 
bias  in  our  predictions,  as  we  are  not  assuming  a  statistical  model  that  we  are  testing  at 
the  same  time.  In  order  to  do  this  however  we  need  to  trust  that  the  wave  tank  accurately 
reproduced  the  wave  spectrum  that  was  assumed  for  the  calculations,  and  that  the 
computational  model  provides  accurate  predictions  of  the  response  on  a  wave  by  wave 
basis.  In  other  words:  We  have  to  assume  that  the  model  uncertainty  added  will  not  affect 
our  results. 


Computational  Model 

The  model  of  the  structural  response  (i.e.  the  horizontal  offset  at  54.8m)  assumes  a  rigid 
body  with  2  degrees  of  freedom  (surge  and  pitch).  The  first-  and  second-order  transfer 
functions  were  established  by  Kim  and  Yue  (1989, 1991).  These  results  were  the  basis  of 
Jha’s  model,  which  includes  additional  effects,  such  as  wave  drift  damping,  viscous 
effects,  mooring  stiffness,  etc.  The  final  results  consisted  of  4  alternative  models  for  the 
spar  buoy. 

1.  Base-case  model  with  diffraction  forces 

2.  Model  with  diffraction  forces  and  wave-drift  damping 

3.  Model  with  diffraction  forces,  wave-drift  damping  and  viscous  force  from 
undisturbed  waves 
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4.  Model  with  diffraction  forces,  wave-drift  damping  and  viscous  force  from 
disturbed  waves 

Jha  concludes  that,  together  with  diffraction  forces,  viscous  effects  should  be  included  to 
accurately  predict  the  mean  offset  of  the  spar.  As  the  simulations  do  not  provide  us  with 
the  disturbed  wave  surface  near  the  structure,  the  third  model  was  selected.  Fortunately 
the  use  of  the  disturbed  wave  surface  showed  no  significant  improvement  of  the  results 
compared  to  the  results  based  on  the  undisturbed  wave  surface  (Jha,  1997). 

For  time  domain  simulation  we  have  to  solve  the  basic  equation  of  motion  in  the  time 
domain. 

Mx{t)+  Cx{t)+  Kx[t)  =  F{t)  4. 1 


Where: 

x(t)  :  horizontal  offset  at  54.8  m 

C 

:  damping  matrix 

M  :  mass  matrix 

K 

:  stiffness  matrix 

F(t)  :  applied  (hydrodynamic)  force 

The  hydrodynamic  force  is  a  function  of  the  random  wave  process  T|(t),  the  geometry  of 
the  structure,  and  the  water  depth. 
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Waves 


The  wave  environment  for  the  test  is  intended  to  represent  the  lOO-year  extreme 
significant  wave  height  sea  state  for  the  Gulf  of  Mexico.  A  Jonswap  spectrum  with  Hs  = 
13.  Im,  Tp  =  14.0s,  and  a  peakedness  factor  of  2.0  was  used  to  represent  these  conditions. 
The  model  tests  however  showed  an  observed  Hs  of  14m,  which  was  used  in  the 
computer  simulations  instead.  From  the  wave  spectrum  complex  Fourier  amplitudes, 
C(tflic)  are  simulated  using: 


4.2 


The  phase  angles  0k,  are  selected  random  from  a  uniform  distribution  between  0  and  27C. 
The  spectral  amplitudes  Sx(tOk)  are  taken  as  exponentially  distributed  random  variables 
with  a  mean  value  given  by  the  wave  spectrum  at  cOk  (Ude,  1996).  These  “noisy” 
simulations  intend  to  reproduce  Gaussian  ensembles  better  than  more  straightforward 
deterministic  methods.  In  addition  “noisy”  simulation  tends  to  predict  the  variance  of  the 
second  order  response  better  (Ude,  1995). 

The  irregular  waves  can  be  written  in  the  time  domain  as: 

rj{t)  =  Re C{q)^  )cxp{ico^t)  4.3 
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Although  TITOP  uses  linear  waves  as  input,  the  transfer  functions  for  the  hydrodynamic 
forcing  also  add  a  second  order  Stokes  component  to  the  waves.  Therefore  effectively  a 
second  order  random  wave  model  was  used. 

Hydrodynamic  Forcing  Model 

Diffraction  analysis  entails  applying  sinusoidal  waves  of  different  frequencies,  (flic, 
selected  from  the  wave  Fourier  spectrum  to  a  structure,  as  described  in  the  previous 
paragraph  to  the  structure.  The  first-order  forces  are  then  found  at  these  incident  wave 
frequencies  as: 

/,  (<)=  Re  )H,  (®.  )exp(>'ffl.>)  4.4 

where  Hi  is  the  first-order  transfer  function  and  the  C(a)j)  coefficients  are  the  complex 
Fourier  amplitudes.  The  second-order  forces  are  found  at  pairs  of  frequencies  as: 

fi  W  =  Re X E 2  )exp[i(ty„  -C0„)t]  4.5 

where  Ha'  is  referred  to  as  the  difference-frequency  quadratic  transfer  function  (QTF). 

Note  that  in  order  to  find  these  transfer  functions  the  spar  was  allowed  to  float  freely 
(Kim  and  Yue  1989, 1991).  Linear  diffraction  analysis  was  used  to  estimate  the 
frequency-dependent  added  mass  and  damping  for  the  spar.  The  second-order  diffraction 
analysis  is  computationally  intensive  and  limited  to  a  few  frequency  pairs  (here  8x8 
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frequency  grid  ranging  from  0.2  to  1.18  rad/s).  A  surface  spline  fitting  scheme  was  used 
to  interpolate  the  sparse  QTF  to  finer  mesh  used  by  TFPOP  (Jha,  1997). 


The  viscous  forces  are  calculated  with  Morison’s  equation,  integrated  from  the  keel  to  the 
free  surface.  The  drag  force  was  based  on  the  absolute  fluid  velocity  with  an  assumed 
drag  coefficient  of  0.6: 


fo  ^)  =  2)  I 


Wheeler  stretching  (Wheeler,  1970)  was  used  to  determine  the  velocity  of  the  water 
particles  above  the  mean  water  level,  for  which  linear  wave  theory  does  not  provide  a 
solution. 


Structural  Model 

The  spar’s  horizontal  offset  in  the  wave  direction  consists  of  2  components,  a 
displacement  in  the  wave  direction,  and  a  rotation  in  the  plane  of  the  waves,  respectively 
its  surge  and  pitch  modes.  The  structure  is  modeled  as  a  rigid  body  with  2  degrees  of 
freedom,  and  the  following  mass,  stiffness,  and  damping  matrices.  (The  first  mode  is 
surge,  while  the  second  is  pitch.) 

4.7 
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where  Mstr  is  the  total  mass  of  the  spar  and  the  water  entrapped  in  the  center  well,  and 
Madd  is  the  (frequency  dependent)  added  mass  matrix,  reflecting  the  effects  of  waves 
radiated  by  the  oscillating  spar.  The  added  mass  was  obtained  from  diffraction  analysis, 
and  is  assumed  to  be  constant  for  low-frequency  modes. 


M 


sTr 


m  -mZcG 
-mZcc 

2.59-10*%  -2.74  10'°%m 

-2.74-10"’%m  3.91-10'^ 


4.8 


2.71-10*% 

-2.60-10'®%m 


-2.60-10"’%m 
3.20- 10'^  %m^ 
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The  symbols  are  explained  in  table  3.1. 

Note  that  Mstr  and  Madd  are  of  the  same  order  of  magnitude  for  this  large-diameter 
structure. 

The  stiffness  matrix,  derived  from  geometry  considerations  is  found  to  be: 


k  -kZ^ 

-  kZf  kZj  +  % 

1.91-10^A^m‘'  -  2.02- 10’ N 

-  2.02- 10’ N  1.60- 10'°  Nm"' 
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where: 


k,  =nR^Hpg{Z,^  -Zcs)~P8R' 


The  eigenfrequencies  associated  with  these  structural  properties  are  331  s  in  surge  and 
69.9  s  in  pitch,  which  correspond  to  the  frequencies  observed  in  the  response  Fourier 
spectrum  and  the  free  decay  tests.  The  eigenmodes  also  confirm  the  model,  as  they  show 
the  correct  geometric  relation  between  surge  and  pitch  (i.e.  1  unit  of  surge  offset  results 
in  zero  pitch  rotation,  and  1  radian  pitch  rotation  results  in  100.6  m  surge  offset). 


What  remains  to  be  modeled  is  the  damping  of  the  structure,  which  is  the  most 
challenging  part.  From  modal  analysis  we  can  construct  the  damping  matrix  C  as  follows: 


0 


0 

2C  (O 

=  p  p 


4.12 


where  cOs,  and  (Op  are  the  natural  periods,  and  and  are  the  damping  ratios  in  surge 
and  pitch  respectively.  O  represents  the  eigenmatrix. 

The  damping  ratios  were  used  to  calibrate  the  model  to  the  observed  response.  This  was 
done  in  an  iterative-procedure  to  match  the-bandwidth  parameterrS  of  the  response 
Fourier  spectrum. 
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A.  =\rs(fW 


4.13 


where: 

A,,  :  the  n*  spectral  moment  S(f)  :  spectral  density  at  frequency  f 

Alternative  methods  and  pitfalls  are  discussed  by  Jha  (1997).  In  addition  we  may  account 
for  wave  drift  damping.  Although  wave  drift  damping  does  not  represent  a  physical  loss 
of  energy,  it  is  also  an  effect  that  is  proportional  to  the  structure’s  velocity.  In  principle  it 
accounts  for  the  structures  relative  velocity  to  the  water.  On  the  hydrodynamic  forcing 
side  of  the  basic  equation  of  motion  there  will  be  a  term  proportional  to  the  velocity  of 
the  water  relative  to  the  strucmre  squared.  This  term  is  linearized,  and  the  velocity  of  the 
structure  is  brought  to  the  other  side  of  the  equation,  where  it  is  grouped  with  the 
damping  terms.  The  calibrated  values  for  ^  will  therefore  depend  on  whether  we  account 
for  wave  drift  damping  or  not.  In  this  case  we  do,  which  results  in  damping  ratios  of  4.0% 
and  0.5%  in  surge  and  pitch  respectively.  As  wave  drift  damping  is  frequency  dependent 
the  elements  of  C  are  not  reported  here.  The  diagonal  terms  for  wave  drift  damping  were 
calculated  using  SWIM  (MIT).  The  off-diagonal  terms  were  estimated  with  Newman’s 
approximation  (Newman,  1974). 
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The  Data 


Using  the  model  of  the  OTRC  spar  36  hours  of  response  data  were  simulated.  Analogous 
to  the  previous  chapter  we  will  analyze  the  data,  the  distributions  of  peaks,  and  the 
distribution  of  the  global  extremes. 


Moments 

The  moments  of  the  process  and  the  peaks  are  reported  in  table  4.1  for  the  total  36  hours. 
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Table  4.1  Moments  of  the  process  and  its  local  peaks  for  1x36  hrs  of  response  data 


Contrary  to  the  1-hour  results,  the  moments  for  the  processes  seem  to  be  very  Gaussian, 
with  skewnesses  near  zero  and  kurtosises  very  close  the  3.  The  peaks  of  each  of  the 
components  also  appear  to  be  well  modeled  by  Rayleigh  distributions,  with  coefficients 
of  variation.close  to  0.52,  and  skewnesses  near  0.63.  Gnly-the-local  peaks  of  the  total 
response,  arising  from  cycles  over  3  separate  time  scales,  depart  clearly  from  the 
Rayleigh  model.  This  may  be  expected,  as  the  total  response  no  longer  shows  the 
narrow-band  frequency  content  consistent  with  the  Rayleigh  model. 
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Besides  the  first-order  wave  frequency  response,  also  surge  and  pitch,  the  second-order 
part  of  the  system’s  response  seem  to  be  linear.  There  are  various  possible  explanations 
for  the  Gaussian  behavior  during  the  simulation,  compared  to  the  model  test,  which  does 
not  seem  Gaussian.  The  computer  model  assumes  the  moorings  to  behave  as  a  linear 
spring,  which  is  a  crude  approximation  of  reality,  especially  for  large  displacements.  The 
moorings  in  the  model  test  are  not  linear,  and  may  cause  a  certain  amount  of  non- 
Gaussianity.  However,  more  importantly,  the  amount  of  data  used  (e.g.  390  peaks  in 
surge)  has  eliminated  most  of  the  noise  and  bias,  which  was  present  in  the  results  for  a 
single  hour  of  data.  Table  4.2  reports  the  means  of  the  36  1-hour  moments  of  he  data  set, 
which  strongly  suggest  that  using  a  limited  1-hour  duration  of  data  biases  the  skewness 
and  kurtosis  of  the  process  and  the  skewness  and  the  coefficient  of  variation  of  the  local 
peaks.Although  the  moments  are  not  as  low  as  the  moments  of  the  single  hour  of 
observed  data,  which  are  reported  in  table  3.3,  the  limited  amount  of  data  explains  a  good 
deal  of  the  bias.  The  higher  moments  of  surge  and  pitch,  which  have  the  least  amount  of 
peaks,  are  affected  most  by  the  data  reduction. 


Process 

mean 

stdv 

m3 

m4 

Surge 

0.00 

3.79 

-0.01 

2.59 

Pitch 

0.00 

3.76 

-0.01 

2.63 

Wave  fr. 

0.00 

2.61 

0.00 

2.96 

Total 

5.11 

6.02 

0.05 

2.87 

Peaks 

mean 

stdv 

m3 

CoV 

Surge 

4.70 

2.23 

0.25 

0.47 

Pitch 

4.76 

2.28 

0.25 

0.48 

Wavefr. 

3.25 

1.71 

0.55 

0.53 

Total 

5.50 

4.34 

0.95 

0.79 

Table  4.2;  Moments  of  the  process  and  its  local  peaks  for  36x1  hr  of  response  data 
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Distribution  of  Peaks 


Distribution  of  Local  Peaks 

For  practical  reasons  we  cannot  display  the  36  1-hour  distributions  per  component. 
Instead  the  models  were  fitted  to  the  entire  36  hours  of  data  (fig  4. 1-4.4)  The  plotted 
distributions  in  fig  4. 1-4.4,  based  on  36  hours  of  data,  clearly  show  that  if  we  restrict 
attention  to  any  single  component  all  of  the  models  behave  quite  well.  There  is  little 
advantage  at  the  component  level  to  generalize  the  Weibull  through  a  quadratic 
transformation.  When  considering  the  total  response  however,  the  quantity  of  ultimate 
interest,  differences  between  the  models  emerge..  The  ‘observed’  distribution  of  the  total 
shows  a  “comer”  in  the  distribution  around  17  m,  which  the  standard  Weibull  fails  to 
capture.  The  quadratic  Weibull  follows  the  data  better  near  the  T  =  1  hr  to  6  hr  return 
period  extremes,  but  it  is  not  clear  how  the  distribution  behaves  even  further  out  in  the 
tail.  The  Hermite  model  also  appear  to  be  able  to  capture  the  tail  behavior  above  the 
“comer”  well  (i.e.  for  mean  return  periods  of  T  =  1  to  6  hrs),  but  fails  to  capture  the  body 
of  the  distribution.  Recall  that  the  Hermite  model  here  transforms  a  narrow-band 
Gaussian  process  (with  Rayleigh  peaks).  Hence  it  does  not  accurately  model  -nor  was  it 
intended  to  model-  the  less  interesting  lower  fractile  peaks  of  a  wide-band  process. 

These  results  clearly  suggest  that  fitting  a  distribution  to  the  total  response  here  should  be 
done  with  care.  In  general  problems  might  arise  when  a  single  limited-parameter  model, 
such  as  the  Weibull,  is  applied  to  cases- where  actual  peaks  arise  from  a  mixture  of 
different  phenomena.  Here,  the  mixture  is  due  to  response  components  with  different  tim- 
scales.  In  previous  applications  to  jack-up  forces  and  responses,  the  mixture  arose  from 


the  combination  of  drag-  and  inertia-induced  effects.  The  Weibull  model  was  found  to 
show  similar,  though  non-conservative  errors  in  this  case  (Torhaug,  1996). 

Distribution  of  10-minute  global  peaks 

The  Gumbel  model,  applied  to  the  10-minute  global  peaks  shows  similar  accuracy  to  the 
best  of  the  preceding  models  (Hermite  and  quadratic  Weibull)  based  on  local  peaks.  In 
particular,  fine  agreement  is  found  through  the  T  =  1  hour  return  period  response,  even 
for  the  total  response.  Comparing  the  “total  response”  plots  in  Figs.  4.4  and  4.6,  all 
models  except  the  Weibull  accurately  predict  the  observed  1-hour  response  of  about  18 
m,  while  the  Weibull  model  returns  a  1-hour  response  of  roughly  20m  (i.e.  about  10% 
error) 

At  or  near  the  T  =  6  hour  mean  return  period,  one  may  view  Figs.  4.5-4.6  as  suggesting 
that  the  Gumbel  model  becomes  somewhat  conservative.  At  the  same  time,  tme  should  be 
wary  about  inferences  based  on  the  last  few  data  points.  (In  comparing  Figs  4. 1-4.4  with 
4.5-4.6,  one  should  note  that  the  first  two  plot  the  offset,  x,  on  log-scale,  thereby  visually 
deemphasising  deviations  in  the  upper-tails.) 


61 


1 


2 


3 


4  5  6  7  8  9  10  12  14  16 1820 

offset  (m) 


Figure  4.1:  Surge  local  peaks  and  fitted  distributions  plotted  on  Weibull- scale 


Figure  4.3:  Wave  frequency  response  local  peaks  and  fitted  distributions  plotted  on 
Weibull-scale 
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Figure  4.4:  Total  response  local  peaks  and  fitted  distributions  plotted  on  WeibulTscale 
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Figure  4.5:  Distribution  of  JO-minute  global  peaks  plotted  on  Gumbel-scale  (36  hours  of  data) 
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The  Extreme  Value  Distribution  (1-hour  global  peaks) 


Results  in  the  preceding  section  treated  the  36-hour  simulation  as  one  long  data  set,  from 
which  one  can  form  “best”  empirical  distributions  to  compare  with  fitted  models.  In  this 
section  we  return  to  considering  the  consequences  of  having  only  one  hour  of  data  (e.g.  a 
single  model  test)  to  estimate  statistics  of  Xmax-  We  now  treat  the  simulations  as  36  1- 
hour  segments,  from  which  we  can  estimate  both  the  bias  and  uncertainty  that  may  arise 
in  a  single,  randomly  selected  hour. 
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Table  4.3  Predictions  of  I -hour  extremes,  and  their  standard  deviation  based  on  36  simulated 
hours,  no  bootstrapping 


Bias 

We  focus  in  table  4.3  on  estimating  |ixmax»  the  mean  value  of  the  hourly  maximum.  The 
“observed”  mean  max  results  are  simply  averages  of  the  36  simulated  hourly  max  values 
Compared  with  these  unbiased  observed  estimates,  the  various  models  show  rather 
consistent  estimates  of  pxmax-  Biases  are  typically  less  than  5%  (table  4.4);  the  only  case 
of  significant  bias  (15%)  occurs  when  applying  the  Weibull  model  to  the  total  set  of 
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response  peaks.  This  is  entirely  consistent  with  the  results  of  the  previous  section.  Note 
also  that  unlike  the  previous  section,  these  results  rely  on  the  assumption  of  independent 
peaks  when  estimating  Fxinax(x)  (see  eq.1.6)  and  hence  Hxmax-  The  agreement  found  here 
suggests  the  adequacy  of  this  assumption  for  extreme  value  estimation. 


Surge 

Pitch 

WsMBfreq. 

Total 

Gunnbei 

1.00 

1.02 

1.00 

1.00 

Weibull 

1.01 

1.09 

1.02 

1.15 

QWeibull 

0.99 

1.05 

0.99 

1.05 

Hemrite 

0.95 

1.02 

1.00 

a99 

Table  4.4  Bias:  the  ratio  between  estimated  and  observed  mean  hourly  extreme. 


Bootstrap  Bias 

The  following  table  shows  the  ratio  between  the  “observed”  standard  deviation  and  the 
average  of  the  36  bootstrap  standard  deviations.  In  order  to  calculate  the  latter,  peaks  of 
each  simulated  hour  were  resampled  2(X)  times  in  the  same  fashion  in  chapter  3,  to 
predict  as  many  mean  1-hour  extremes.  From  these  2(X)  predictions  the  “bootstrap” 
standard  deviation  is  calculated.  Each  of  the  36  simulated  hours  was  treated  in  this  way  to 
estimate  the  standard  deviation;  the  bootstrap  result  reported  here  is  the  average  of  these 
36  values.  If  bootstrapping  provides  an  unbiased  estimate  of  the  standard  deviation,  this 
ratio  should  be  equal  to  1. Table  4.5  however  shows  this  ratio  is  2.1  on  average,  which 
clearly  suggests  that  bootstrapping  may  underestimate  the  standard  deviation. 
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1.95 

1.55 
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1.93 
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123 

1.79 

QWabiJI 
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243 

1.32 

1.68 

Table  4.5  Ratio  between  the  “observed”  and  bootstrap-estimated  standard  deviation 

Note  that  for  the  Hermite  no  bootstrapping  could  be  performed  as  it  is  fitted  to 
the  process,  while  we  can  only  resample  peaks. 
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Table  4.5  indicates  that  the  bootstrap  estimate  of  the  standard  deviation  is  much  less 
biased  for  the  wave  frequency  response  predictions  made  with  both  Weibull  models.  The 
relatively  large  number  of  wave  frequency  peaks  includes  data  far  in  the  upper  tail,  which 
reduces  the  effect  of  the  truncation  of  the  distribution  from  which  synthetic  data  sets  are 
sampled.  This  truncation  effect,  which  we  believe  may  bias  our  estimate  of  the  standard 
deviation,  is  much  more  severe  for  surge  and  pitch  as  they  have  less  peaks.  Surge  and 
pitch  peaks  will  also  show  higher  correlation  between  peaks,  an  effect  which 
bootstrapping  analysis  ignores. 


Uncertainty 

The  wave  frequency  response  also  shows  some  reduction  of  the  standard  deviation  of  our 
prediction  using  a  model  (table  4.6).  This  effect  is  not  observed  for  pitch  and  surge. 

Again  we  expect  this  reduction  to  be  the  result  of  the  relatively  large  number  of  peaks  in 
our  data  set,  which  reduces  uncertainty  in  our  estimate  of  the  distribution  parameters. 
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Harrite 

0.91 

1.13 
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Table  4.6  Ratio  between  the  “observed”  and  the  from  a  model  predicted  standard 
deviation  of  the  I -hour  mean  extreme  response 


An  interesting  observation  is  also  that  coefficient  of  variation  of  the  peaks  of  the  second 
order  components,  surge  and  pitch,  is  approximately  twice  as  large  as  the  coefficient  of 
variation  of  the  peaks  of  the  first  order  component,  the  wave  frequency  response.  This  is 
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an  effect  that  is  generally  observed  for  first  and  second  order  processes  (e.g.  Torhaug 
1996). 


SRSS 

As  it  appears  from  the  previous  paragraph  that  the  distribution  of  the  total  response  is 
difficult  to  model,  one  might  wonder  whether  the  mean  extreme  total  value  can  be 
estimated  from  its  components,  for  which  fairly  straightforward  models,  such  as  a 
Rayleigh  distribution,  seem  to  suffice.  The  standard  method  that  is  used  in  random 
vibrations  for  this  is  the  “square  root  of  the  sum  of  squares  method”,  or  srss.  Although  the 
method  does  not  have  a  strong  analytical  background,  it  often  gives  good  results  for 
uncorrelated  Gaussian  processes.  Extremes  of  a  Gaussian  process  tend  to  be  proportional 
to  the  standard  deviation  of  the  process.  The  standard  deviation  of  the  non-Gaussian  total 
is  equal  to  the  “square  root  of  the  sum  of  squares”  of  the  standard  deviations  of  the 
Gaussian  (Rayleigh)  components,  (plus  the  correlation  terms  which  may  become 
important  for  closely  spaced  modes).  The  prediction  of  the  mean  total  1-hour  extreme 
response  using  the  srss  method  for  the  observed  values  is  V(8.76^  +  8.96^  +  9.83^)  = 

15.93  m.  The  observed  value  is  19.12  m.  This  relatively  large  difference  may  be  due  to 
either  the  correlation  among  these  components,  or  to  their  non-Gaussian  behavior  (Naess 
and  R0yset  1998).  The  srss  method  does  not  seem  very  useful  in  this  particular  case. 
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Improving  Use  of  Limited  Data 


In  the  previous  sections  we  concluded  that  fitting  a  model  does  not  always  reduce  the 
uncertainty  in  our  predictions  of  the  mean  (and  distribution)  of  the  hourly  extreme  value. 
All  of  our  analysis  however  has  been  concentrated  on  predictions  for,  and  based  on,  a  1- 
hour  time  period.  Previous  work  by  Torhaug  (1996)  shows  that  uncertainty  can  be 
reduced  by  a  factor  of  roughly  2.5  for  a  jack-up  structure,  by  fitting  a  model.  In  order  to 
make  a  comparison  between  the  statistics  of  the  extreme  response  of  a  spar  buoy  and  a 
jack-up  structure,  we  will  mimic  the  results  found  by  Torhaug  in  this  section.  The  most 
important  difference  is  that  the  jack-up  results  were  based  on  using  1-hour  response 
histories  to  predict  6-hour  mean  extreme  responses.  The  following  table  shows  the  results 
for  the  6  observed  6-hour  extremes,  and  the  36  predicted  mean  6-hour  extremes,  based  on 
1  hour. 
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Surge 


Wave  frequency  response 
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Table  4.7  Predictions  of  6-hour  extremes,  and  their  standard  deviation  based  on  36  simulated  hours,  no 
bootstrapping 


Bias 

The  ratio  between  the  best  estimate  of  the  mean  6-hour  extreme  value  and  the  observed 
6-hour  extreme  is  reported  in  table  4.8.  The  differences  between  the  predicted  and 
observed  values  are  often  more  pronounced  because  we  are  extrapolating  a  given  model 
further  into  its  tail.  Now  the  Hermite  predictions  seem  relatively  low,  while  the  Gumbel 
predictions  appear  to  remain  nearly  unbiased.  This  may  be  because  the  Gumbel  model  is 
already  fit  to  10-minute  maxes,  hence  requires  the  least  extrapolation  to  reach  the  6-hour 
level. 
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1.02 

1.01 

1.00 

Weibull 
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Table  4.8  Ratio  of  predicted  extremes  to  observed  6-hour  extreme 
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Uncertainty 

Table  4.7  shows  the  standard  deviation  increases  if  we  estimate  the  mean  6-hour  extreme 
not  with  the  observed  6-hour  extreme,  but  rather  with  a  model  fit  prediction  that  only 
uses  1  hour  of  data.  This  effect  is  not  surprising:  results  based  on  1  hour  of  data  should  be 
more  variable  than  those  based  on  6  hours.  What  is  interesting  however  is  that  the 
increase  in  standard  deviation  of  a  single  prediction  of  the  6-hour  extreme  is  only  roughly 
a  factor  of  1.5  on  average. 

What  is  reported  in  the  first  column  of  table  4.7  is  the  mean  of  36  of  these  estimates  for 
the  fitted  models,  and  the  mean  of  6  observed  6-hour  extremes.  The  standard  deviation  of 
these  estimated  means,  due  to  having  a  limited  sample  of  36  hours  total,  is  the  standard 
deviation  reported  in  table  4.7  divided  by  the  square  root  of  the  number  of  estimates.  As 
the  number  of  estimates  is  6  for  the  observed  extremes  and  36  for  the  predicted  extremes, 
the  ratio  between  the  observed  and  fitted  standard  deviations  decreases  by  the  square  root 
of  6.  The  standard  deviation  of  our  best  estimate  of  the  6-hour  extreme  therefore 
decreases  by  a  factor  V6/1.5  =  1.6  on  average,  suggesting  extrapolating  does  help, 
although  it  is  not  as  high  as  the  2.5  observed  by  Torhaug.  Apparently  the  increase  of 
uncertainty  by  extrapolating  is  less  than  the  reduction  caused  by  using  more  data. 

The  reductions  of  the  coefficient  of  variation  are  reported  below. 
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Table  4.9  Ratio  of  the  model-predicted  CoV  of  the  mean  6-hour  extremes  to  the  CoV  of 
the  6  observed  extremes 


Bootstrap  Bias 

The  following  table  confirms  the  bias  of  the  bootstrap  estimate  of  the  standard  deviation 
again,  which  is  higher  as  we  are  extrapolating  further  out  into  the  tails  of  the  underlying 
distributions.  Moreover  the  truncation  of  the  distribution  of  peaks  now  also  affects  the 
predictions  for  the  wave  frequency  response. 
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Table  4.10  Ratio  of  observed  CoV  to  bootstrap-predicted  CoV  of  6-hour  extremes 


Extrapolation 

The  previous  results  for  extrapolating  leave  us  with  the  question  of  how  small  should  we 
make  the  time  periods  to  which  we  fit  a  distribution  of  extremes  in  order  to  predict  an 
extreme  value  with  the  least  amount  of  uncertainty.  The  smaller  the  period,  the  more  we 
leverage  the  effect  of  having  more  estimates.  On  the  other  hand  the  shorter  the  periods, 
the  further  we  need  to  extrapolate.  The  following  tables  report  the  bias  and  the  reduction 
of  the  coefficient  of  variation  compared  to  1  hour  for  fitting  the  models  to  20-minute 
periods,  producing  108  estimates  of  the  mean  hourly  and  6-hour  extreme. 
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Table  4.11  Bias  in  1-hour  mean  max  predictions,  based  on  20  mins 
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Table  4.12  Bias  in  6-hour  mean  max  predictions,  based  on  20  mins 


Not  surprisingly  relatively  small  biases  become  fairly  large  when  we  extrapolate  further 
out  in  the  tail.  As  the  limited  amount  of  data  to  which  the  models  are  fit  can  bias  the 
higher  moments  relatively  low,  these  higher-moment  models  may  tend  to  underpredict 
extreme  responses. 


The  effect  of  this  extrapolation  on  the  uncertainty  in  our  estimates  is  shown  in  tables  4.13 
and  4.14.  These  report  the  ratio  of  the  coefficient  of  variation  of  the  108  1-hour  and  6- 
hour  predictions  of  the  mean  max  based  on  20  minutes,  to  the  mean  “observed”  standard 
deviations  for  1-  and  6-hour  periods  (symbolically:  Oes/CTdirect.  for  1  and  6  hour  mean 
extremes).  Table  4.13  is  directly  comparable  to  table  4.6. 
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Table  4.13  CoV  reduction  of  I -hour  mean  max  based  on  20  mins  vs.  I  hr 
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Table  4.14 

i  CoV  reduction  of6’hour  mean 

max  based  on  20  mins  vs,  6hrs 
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Without  adding  much  extra  bias,  the  variability  in  our  estimate  of  the  mean  reduces 
despite  extrapolating  from  a  small  amount  of  data  for  surge  in  particular  (<  4  values  on 
average).  The  average  CoV  reduction  for  the  components  is  0.76  (a  factor  of  1.3).  We 
cannot  reduce  the  time  segment  length  any  further  for  surge  in  view  of  its  natural  period 
of  5  minutes,  as  we  need  at  least  3  peaks  per  period  to  fit  the  models.  A  similar  argument 
has  to  be  made  for  the  Gumbel  distribution,  which  needs  at  least  2  global  periods  within 
the  time  segment  that  will  be  fitted  (20  mins).  In  a  global  period  of  10  mins  there  will  be 
2  peaks  in  surge  on  average,  for  which  a  local  peak  model  would  probably  be  more 
appropriate.  A  global  period  less  than  5  minutes  would  be  meaningless  for  surge.  Based 
on  the  previous  considerations  the  segment  length  of  20  minutes  was  regarded  as  an 
absolute  minimum  for  surge  in  this  case.  For  the  other  components  the  data  are  not  as 
scarce.  The  following  tables  show  the  result  of  extrapolating  6-minute  periods  to  1-hour 
extremes,  and  can  be  compared  directly  to  tables  4.1 1  and  4.13. 


Suge 

Rtch 

Wawefr. 

Total 

Gurtel 

- 

0.97 

1,00 

1.22 

WeibJI 

- 

0.06 

1.01 

1.15 

QWeitxJI 

- 

0.66 

078 

071 

Hermite 

- 

077 

094 

082 

Table  4.15  Bias  in  1-hour  mean  max  predictions,  based  on  6  mins 


Suge 

Pitch 

Wawefr. 

Total 

GUrtel 

062 

069 

075 

WeibLil 

- 

0.64 

0.50 

0.66 

QWeibiil 

- 

062 

0.51 

054 

Harrrite 

- 

063 

063 

054 

Table  4.16  CoV  reduction  of  1-hour  mean  max  based  on  6  mins  vs.  Ihr 
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As  expected  the  data  are  insufficient  for  l-hour  extreme  predictions  in  surge.  However, 
also  for  pitch  and  the  total  and  to  a  lesser  extent  for  the  wave  frequency  response,  the 
results  are  very  biased  compared  to  the  previous  results.  There  still  is  a  reduction  of  the 
CoV  of  our  predictions,  but  this  irrelevant  with  the  amount  of  bias  observed.  Fitting  to  6- 
minute  periods  and  extrapolating  to  10-fold  periods  of  1  hour  clearly  shows  that  using 
periods  that  are  too  short  gives  very  poor  results. 


Conclusions 

The  results  suggest  that  it  is  possible  to  reduce  the  uncertainty  in  our  estimates  of  the 
mean  extreme  response,  from  limited  data,  by  producing  as  many  estimates  of  the  mean 
extreme  response  as  possible  by  fitting  a  model  to  short  time  segments.  The  limiting 
factors  for  the  duration  of  these  time  segments  is  the  minimum  amount  of  data  required  to 
estimate  the  moments  to  fit  to,  and  the  amount  of  bias  introduced  in  our  predictions  by 
using  less  data. 

The  results  also  suggest  that  the  moments  for  the  single  hour  of  model  test  data  are  biased 
due  to  limited  data,  when  they  are  compared  to  data  from  36  hours  of  computer 
simulation.  The  local  peaks  of  the  components  of  the  36  hours  of  data  showed  very  good 
agreement  with  a  Rayleigh  model.  The  local  peaks  of  the  total  response,  which  is  the  sum 
of  the  three  components,  seems  to  be  better  described  by  a  Hermite  or  quadratic  Weibull 
model.  The  Gumbel  model  of  the  10-minute  global  peaks  showed  good  agreement  for  the 
total  and  all  components. 


78 


1  Conclusions 


In  this  report  we  have  discussed  estimation  methods  of  the  statistics  of  the  extreme 
response  of  a  spar  buoy,  during  a  lOO-year  Gulf  of  Mexico  extreme  sea  state.  We 
analyzed  the  statistics  of  a  1-hour  time  series  of  the  horizontal  offset  at  54.8  m  MWL, 
from  a  model  test,  and  simulated  36  hours  of  response  data  using  a  computer  model 
calibrated  against  the  same  model  test.  The  statistics  of  interest  were  the  mean  of  the 
distribution  of  the  extreme  response  during  a  certain  period  of  time,  and  the  standard 
deviation  of  our  estimate  of  that  mean,  quantifying  the  accuracy  of  our  estimate.  The 
results  were  compared  and  suggestions  are  made  for  more  effective  use  of  limited  data. 

The  total  spar  buoy  response  consisted  of  3  narrow-band  processes  with  dominant 
frequencies  at  the  structures  natural  periods  in  surge  and  pitch,  and  at  the  spectral  peak 
period  of  the  sea  state.  Predictions  of  the  mean  extreme  hourly  and  6-hourly  response 
were  made  for  the  total  response  zind  all  3  components. 

Predictions  for  a  particular  period  of  time  were  made  fitting  a  parametric  distribution  to 
that  period,  using  the  10-minute  global  peaks,  all  the  local  peaks,  and  the  entire  process, 
in  order  to  investigate  the  effect  of  the  amount  of  data  included  in  the  fit.  The 
distributions  were  fitted  to  the  moments  of  the  selected  data.  The  distribution  models 
used  were:  (1)  the  Gumbel  model  (fitted  to  the  l"  2  moments  of  the  10-minute  global 
peaks),  (2)  the  standard  Weibull  model  (fitted  to  1*‘  2  moments  of  the  local  peaks),  (3)  a 
quadratic  Weibull  model  (fitted  to  the  1*‘  3  moments  of  the  local  peaks),  and  (4)  a 
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Heraiite  transformation  of  the  Rayleigh  model  (fitted  to  the  1*'  4  moments  of  the 
process). 

From  the  long  duration  simulated  data  (36  hours)  we  concluded  that  each  of  the  3 
frequency  components,  the  2-parameter  models  capture  their  apparently  Gaussian 
narrow-band  behavior.  For  the  observed  short  duration  data  (1  hour)  the  distributions  of 
the  components  seemed  considerably  more  narrow  and  “pinched”.  We  believe  this  to  be 
largely  due  to  the  effect  of  limited  data,  which  biases  the  higher  moments  in  particular.  A 
similar  bias  was  found  from  our  36  hours  of  simulation,  if  we  processed  it  to  estimate 
moments  separately  for  each  hour.  Of  course  limitations  of  the  computational  model  may 
also  account  partly  for  this  discrepancy.  The  model  used  here  was  based  on  2"“*  order 
diffraction  analysis,  including  wave  drift  damping,  viscous  forces,  and  a  constant 
mooring  stiffness.  The  total  response  is  wide-band  mixture  of  3  component  processes, 
where  peaks  are  not  captured  well  by  the  standard  Weibull  models.  This  causes  bias  in 
our  predictions  of  the  statistics  of  the  extreme  response.  The  Hermite  model,  quadratic 
Weibull  model,  and  the  Gumbel  model  matched  the  total  response  better,  and  showed 
very  little  bias. 

Estimates  of  the  uncertainty  in  our  predictions  were  made  using  non-parametric 
bootstrapping,  and  alternatively  through  simulation  (or  parametric  bootstrapping).  The 
results  showed  that  bootstrapping  local  peaks,  and  negelecting  their  correlation,  may  not 
be  very  suitable  to  estimate  statistics  of  extremes.  As  compared  with  the  uncertainty 
found  by  long  simulations,  bootstrapping  was  found  to  underestimate  the  standard 
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deviation  of  our  predictions  of  the  mean  1-hour  and  6-hour  extremes.  We  therefore 
encourage  the  use  of  simulation  -  as  opposed  to  non-parametric  bootstrapping  -  to 
consider  these  issues  of  bias  and  uncertainty.  While  the  computational  model  may  not  be 
perfect,  it  serves  as  a  useful  basis  to  preserve  the  distribution  and  correlation  structure  of 
the  data,  based  on  the  best  physical  information  at  hand. 

Based  on  the  results  of  36  hours  of  simulated  response  it  appears  that  for  all  models  the 
uncertainty  in  the  predictions  of  the  mean  extreme  value  can  be  reduced  by  dividing  the 
limited  data  in  short-duration  segments  from  which  many  estimates  of  the  mean  extreme 
response  can  be  made  by  extrapolating  a  fitted  model.  It  appeared  that  relatively  very 
short  segments  reduced  uncertainty  most,  without  adding  much  bias. 
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Appendix  A: 


Motions  of  a  Spar  Buoy  in  Random  Seas: 
Comparing  Predictions  and  Model  Tests 

Alok  K.  Jha,  P.R.  de  Jong,  and  Steven  R.  Winterstein 
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Motions  of  a  Spar  Buoy  in  Random  Seas: 
Comparing  Predictions  and  Model  Test  Results 

Alok  K.  Jha,  R  R.  de  Jong,  and  Steven  R.  Winterstein 


Civil  Engineering  Dept.,  Stanford  University 


Abstract 

This  study  compares  the  analytically  predicted  motions  of  a  floating  spar  buoy  platform  with  the 
results  of  wave  tank  experiments.  Results  studied  include  extreme  conditions  in  both  the  Gulf  of 
Mexico  and  the  North  Sea.  Base-case  predictions  combine  nonlinear  diffraction  loads  and  a  linear, 
multi-degree-of-freedom  model  of  the  spar  stiffness  and  damping  characteristics.  Refined  models 
add  the  effect  of  wave-drift  damping,  and  of  viscous  forces  as  well.  Consistent  choices  of  damping 
and  wave  input  are  considered  in  some  detail.  These  successive  model  refinements  are  generally 
found  to  improve  agreement  with  the  model  test  results. 

Keywords 


Nonlinear  wave  diflFraction;  offshore  structures;  random  vibration;  spar  buoys,  structural  reliability; 
wave  tank  experiments. 

Introduction 


This  study  describes  ongoing  research  into  the  statistical  response  behavior  and  reliability  of  a 
particular  deep-water  floating  structure:  a  spar  buoy.  Typically,  the  spar  buoy  concept  involves  a 
deep-draft,  large-diameter  cylindrical  floating  structure,  with  slack  or  taut  mooring  (e.g.,  Glanville 
et  al,  1991).  This  concept  has  recently  gained  increased  interest  within  the  offshore  community.  For 
example.  Oryx  has  installed  the  first  production  spar  (Neptune)  in  1996,  while  Chevron  is  currently 
designing  the  first  spar  (Genesis)  for  both  drilling  and  production.  Concurrently,  a  particular 
spar  buoy  has  been  designated  the  “theme  structure”  of  the  NSF-sponsored  Offshore  Technology 
Research  Center  (OTRC),  centered  at  Texas  A&M  University  and  at  the  University  of  Texzis  at 
Austin. 

This  study  compares  analytical  predictions  of  spar  surge  motions  with  the  results  of  model  tests  from 
the  OTRC  wave  tank  (OTRC,  1995).  Comparisons  are  shown  both  for  summary  response  statistics 
and  for  complete  time  histories.  Consistent  choices  of  damping  and  wave  input  are  considered  in 
some  detail.  Responses  are  filtered  and  compared  for  three  distinct  frequency  ranges:  a  relatively 
high-frequency  contribution  due  to  first-order  wave  energj^,  a  low-frequency  contribution  due  to 
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Figure  1:  Elevation  view  of  spar  platform. 
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Figure  2:  Degrees  of  freedom  for  spar. 


pitch,  and  a  still  lower  frequency  contribution  due  to  surge.  Model  tests  are  studied  for  extreme 
(100-year)  conditions  in  both  the  Gulf  of  Mexico  and  the  North  Sea. 

Results  from  analytical  models  are  shown  here  over  a  range  of  increasing  modelling  detail.  The  base 
case  includes  nonlinear  diffraction  forces  and  a  linear,  multi-degree-of-freedom  structural  model. 
Refinements  on  this  include  the  addition  of  wave  drift  damping,  and  then  of  viscous  forces  as  well. 
These  successively  more  detailed  models  are  generally  found  to  yield  improved  agreement  with 
model  test  results.  The  analytical  predictions  also  show  the  ability  to  capture  another  notable 
feature  of  the  spar  model  tests;  namely,  the  apparent  “mode-swapping,”  between  the  spar  response 
in  pitch  and  surge  modes,  during  the  hour-long  tests. 

Spar  Buoy  Characteristics 

Figure  1  shows  the  prototype  dimensions  of  the  spar  buoy  under  study.  Note  its  relatively  deep 
draft  (i/=198.2m),  particularly  with  respect  to  its  diameter  (D=40.5m).  For  prediction  purposes 
the  spar  buoy  hull  is  assumed  rigid,  and  its  mooring  lines  are  modelled  as  a  set  of  massless,  linear 
springs.  To  predict  the  spar’s  motions  in  the  along-wave  direction,  we  adopt  a  2DOF  model  that 
includes  the  surge  motion  Xi(i)  and  pitch  rotation  X5(t)  at  the  mean  water  level*  (Figure  2).  At  an 
elevation  z  above  this  level,  the  corresponding  along-wave  displacement  of  the  rigid  spar  is  predicted 
simply  as  Xi (t)  4“  2:  •  x^{t).  In  particular,  we  apply  this  result  here  with  z=54.8m,  to  compare  with 
video- recorded  surge  motions  at  this  elevation  during  the  wave  tests  (OTRC,  1995). 

Mode  Shapes  and  Periods 

Assuming  small  deformations,  the  2x2  stiffness  and  mass  matrices  can  be  constructed  from  geomet¬ 
rical  considerations  (Jha,  1997).  This  mass  matrix  includes  added  mass  terms,  reflecting  first-order 

*Note  that  we  retain -the  connnon  convention  that  numbers^urge  and  pitch  DOFs  as '**1”  and  “5”  respectively, 
although  no  other  DOFs  are  included  here. 


84 


(1) 


wave  radiation  effects.  The  resulting  mode  shapes  and  natural  frequencies  are 
/i  =  ^[HzUi  =  [l  0]^;  /s  =  1  [HzUs  =  (100  if 

These  modal  frequencies  agree  well  with  the  natural  periods,  Ti=330s  and  T5=67s,  estimated  from 
free-decay  tests  of  the  spar  (OTRC,  1995).  Note  that  this  lower-frequency  mode  involves  a  pure 
translation,  while  the  higher-frequency  mode  reflects  a  pure  rotation  about  an  axis  located  at  depth 
100m  below  the  MWL.  (Equivalently,  Eq.  1  implies  that  a  small  rotation  xs  [rad]  is  accompanied 
by  a  translation  of  xi=100x5  [m]  at  the  MWL.)  These  modes  directly  reflect  the  translational  and 
rotational  stiffnesses,  respectively,  of  the  spar’s  mooring  system. 


Qualitative  Results  and  Consistent  Damping  Estimates 

We  consider  here  the  spar  model  tests  that  reflect  extreme,  roughly  100-year  wave  conditions.  We 
also  focus  on  tests  that  apply  wave  loads  onlj^  neglecting  other  tests  that  include  simultaneous 
current  and/or  wind  loads.  This  leaves  us  with  three  model  tests,  each  lasting  1  hour  (all  time  and 
length  units  here  reflect  prototype  scale).  Two  of  the  three  are  separate  realizations  of  100-year 
Gulf  of  Mexico  seastates,  while  the  third  models  100-year  North  Sea  conditions.  We  refer  here 
to  these  seastates  as  “GOMl”,  “GOM2”,  and  “NS”.  (In  OTRC  internal  reporting,  these  tests  are 
respectively  denoted  “aran3”,  “aran4”,  and  “aran5”.  To  date,  only  “aran3”  has  received  systematic 
study  by  OTRC  investigators;  e.g..  Ran  et  al,  1996,  Weggel  and  Roesset,  1996). 

Wave  Measurements  and  Characteristics 

A  reference,  “undisturbed”  wave  elevation  history  has  been  measured  during  the  OTRC  tests  by 
a  probe  located  125m  (prototype  scale)  from  the  spar,  in  a  direction  perpendicular  to  the  wave 
direction.  The  spectra  of  these  waves  are  found  to  be  relatively  well-fit  by  JONSWAP  spectral 
shapes  with  7=2;  the  significant  wave  height  and  peak  period  values  are  estimated  as  /75=14.1m 
and  Tp=14.1s  for  the  Gulf  of  Mexico  seastates,  and  Hs=l^-8m  and  7^=16.15  for  the  North  Sea  test 
(Jha,  1997).  Note  however  that  our  response  predictions  use  the  observed  wave  histories  from  the 
tests,  and  not  the  simulated  input  from  a  theoretical  wave  spectral  model. 

Response  Measurements  and  Characteristics 

Figure  3  shows  the  power  spectrum  of  the  spar  displacement,  measured  at  height  2=54.8m  above 
MWL,  during  the  GOMl  test.  Note  its  two  low-frequency  modes,  at  around  /i=l/330  and  /5=l/70 
Hz,  reflecting  motions  induced  by  surge  and  pitch  resonance.  As  Figure  3  shows,  we  use  bandpass 
filters  here  to  separate  the  observed  surge  component  (0-.006  Hz),  pitch  component  (.006”.03  Hz), 
and  remaining  wave  frequency  component  (above  .03  Hz).  This  gives  rms  response  contributions 
of  <T5urpe=3.4m,  apiicA=4.0m,  and  au;ave=2.5m.  Thus  a  linear  force  model,  which  predicts  energy 
only  at  the  wave  frequencies,  would  capture  only  a  small  portion  of  the  response  rms.  It  would  also 
completely  fail  to  predict  the  mean  response,  here  found  to  be  4.9m.  The  other  tests  offer  similar 
results.  This  shows  the  need  for  models  of  nonlinear  forces — diffraction,  drag  or  both — to  explain 
not  only  the  mean  offset  but  also  the  amplitude  of  slow-drift  oscillations  for  the  spar.  The  effects 
of  both  nonlinear  diffraction  and  drag  loads  are  considered  below. 

Figure  4  shows  the  corresponding  time  history  of  the  response  during  the  GOMl  test.  Both  the  total 
response  and  its  filtered  components  are  shown.  Consistent  with  its  power  spectrum  in  Figure  3, 
the  response  indeed  displays  three  distinct  time  scales.  What  Figure  3  fails  to  reveal,  however,  is 
that  the  relative  contribution  of  the  different  frequency  components  does  not  remain  constant  over 
time.  For  the  GOMl  test,  the  observed  response  changes  qualitatively  at  around  t=1500s,  when 
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GOM1 :  Total  Surge  Spectrum 


Frequency  (Hz) 


Figure  3:  Spectrum  of  measured  surge  displacements  at  54.8m  elevation  above  MWL  for  GOMl 
seastate 


the  surge  component  begins  to  contribute  significantly.  The  other  tests  show  somewhat  similar 
shifts  between  the  energ>^  in  surge  and  pitch  modes — although  this  “mode  swapping”  is  observed 
at  different  times,  and  for  different  durations,  in  different  tests.  The  wave  input  histories  show  no 
such  episodic  nature.  This  reflects  a  further  modelling  challenge:  can  analytical  models  predict 
not  only  the  correct  average  frequency  content  (e.g.,  the  spectrum  in  Figure  3),  but  also  time- 
domain  behavior  consistent  with  Figure  4?  Clearly,  this  time-domain  evolution  of  surge  and  pitch 
components  depends  directly  on  (1)  their  initial  conditions  at  the  beginning  of  the  test  recording 
and  (2)  the  damping  values  assigned  to  these  modes.  We  therefore  discuss  these  issues,  particularly 
damping  estimation,  in  some  detail. 

Estimating  Initial  Conditions 

In  the  experiments,  spar  motions  were  recorded  after  about  15  minutes  (prototype  scale),  when  the 
wave  tank  conditions  were  deemed  to  have  achieved  steady-state  conditions.  Thus,  the  assumption 
of  at-rest  initial  conditions  would  corrupt  our  predictions,  more  so  in  the  surge  mode  which  contains 
relatively  few  cycles  over  the  hour-long  test.  To  avoid  this,  our  predictions  use  initial  conditions 
consistent  with  the  tests;  i.e.,  for  each  test  we  filter  the  observed  motions  to  estimate  surge  and 
pitch  components  (e.g..  Figure  4).  The  initial  values/velocities  of  these  components  are  then  used 
to  start  our  slow-drift  motion  predictions  (Jha,  1997). 

Estimating  Surge  and  Pitch  Damping 

Because  the  tests. include  relatively  few -cycles  of  lightly  damped  motion,  it  is  challenging  to  form 
precise  damping  estimates  from  them.  We  focus  here  on  frequency-domain  damping  estimates,  using 
response  spectra  from  the  various  tests  (e.g..  Figure  3).  One  may,  for  example,  select  dampings 
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Time  (sec) 

Figure  4:  Total  measured  horizontal  displacement  and  its  filtered  surge,  pitch  and  wave-frequency 
components  for  GOMl  seastate 


Q  so  that  our  analytical  model  predicts  the  correct  area  under  each  of  these  observed  spectral 
modes  (i.e.,  the  variances  cr'^^rge  ^pitch)‘  ^  danger  in  this  approach  is  that  it  may  mask  a 
force  modelling  error  (e.g.,  Ude,  1994).  For  example,  if  predicted  forces  are  overestimated  we 
may  overestimate  damping  in  an  effort  to  compensate.  Thus,  we  instead  seek  damping  values 
to  preserve  the  observed  spectral  shape — ^the  rms  values  <j surge  and  Cpuch  are  reserved  to  test  the 
model’s  adequacy.  One  such  measure  of  spectral  shape  is  the  half-power  bandwidth,  at  which  the 
response  spectrum  decays  to  half  its  peak  value.  By  setting  this  bandwidth  to  its  approximate  value 
from  theory — fnp—^Qfi  to  either  side  of  the  natural  frequency  /,* — one  can  estimate  the  damping 
Q  from  an  obser\^ed  power  spectrum.  For  example,  if  the  spar  buoy  has  damping  Ci=-05  at  the 
surge  frequency  /i= 1/330,  we  find  /i/p=l/6600  Hz,  Unfort unatel}^  from  a  T=l-hour  history,  our 
finest  frequency  resolution  is  d/=l/r— 1/3600  Hz — too  coarse  to  resolve  the  half-power  bandwidth, 
even  if  no  frequency-averaging  is  applied  to  the  observed  spectrum. 

We  are  therefore  led  to  consider  the  average  shape  of  the  response  spectrum  across  frequencies,  as 
measured  by  the  unitless  bandwidth  measure  6  (Vanmarcke,  1972): 

<5  =  ^1  -  A?/(AoA2)  ;  =  I  rs,if)df  (2) 

Note  that  in  general,  any  parameter  of  the  form  6n=[l  -  could  be  used  to  reflect 

bandwidth:  -4  0  as  the  bandwidth  narrows.  Perhaps  the  most  widely  used  is  62 ^  e.g,  in  modelling 

peaks  of  a  Gaussian  process.  We  use  S  with  n=l  here,  as  its  lower  spectral  moments  are  less  sensitive 
to  high-frequency  spectral  content.  W^e  apply  Eq.  2  twice,  over  the  frequency  ranges  of  surge  (0-.006 
Hz)  and  pitch  (.006-.03  Hz)  components,  to  find  separate  S  values  that  characterize  their  respective 
modal  bandwidths. 

For  a  IDOF  system  under  broad-band  loads,  6  can  be  related  directly  to  the  damping  level  (Van¬ 
marcke,  1972,  Ude  and  W^interstein,  1996).  To  form  estimates  consistent  with  our  2DOF  model, 
however,  we  select  damping  ratios  Ci  and  Q  so  that  our  predicted  response  shows  the  same  S  val¬ 
ues,  in  both  the  surge  and  pitch  frequency  ranges,  as  found  from  the  observed  responses.  This  is 
an  iterative  process,  which  must  be  performed  for  each  choice  of  (1)  seastate  and  (2)  predictive 
response  model.  We  differentiate  here  between  2  seastates  (GOMl  and  GOM2  versus  NS),  and 
among  4  predictive  models.  These  models  are  described  below,  together  with  the  2x4  values  of 
modal  dampings  that  result.  In  each  case,  a  corresponding  damping  matrix  C  is  inferred  from  the 
mass  matrix  and  the  matrix  $  of  modal  shapes:  in  terms  of  (?=diag(47rG/i). 

Predictive  Models  of  Forces  and  Damping 

Our  first,  base-case  model  applies  diffraction  forces  only  (the  “DF”  model).  Linear  diffraction 
gives  first-order  transfer  functions  Fi^^(a;jt)  and  F5^^(a;jt),  the  (complex)  amplitude  of  surge  force 
and  pitch  moment  due  to  a  unit-amplitude  wave  at  frequency  A  corresponding  second-order 
diffraction  analysis  gives  and  the  surge  and  pitch  excitation  amplitudes 

at  the  difference  frequency  —  ujj  due  to  pairs  of  incident  waves  at  frequencies  cVfc  and  cvj  (Kim 
and  Yue,  1989;  Kim  and  Yue,  1991),  These  subharmonic  excitations  drive  the  surge  and  pitch 
resonant  motions,  which  dominate  the  spar  response  (e.g.,  Figure  3).  Note  that  these  diflFraction 
forces  assume  as  input  not  the  undisturbed  total  wave  elevation/potential  but  rather  its  first-order 
component.  Here  we  use  new  methods  (W^interstein  and  Jha,  1997)  to  identify  the  underlying 
first-order  contribution  to  the  observed  undisturbed  wave. 

For  this  model,  the  damping  matrix  C  gives  the  major  source  of  damping.  For  example,  Table  1 
shows  that  this  DF  model  requires  the  damping  ratios  Ci=4.5%  and  ^5=1  *6%  to  match  the  spectral 
bandwidths  estimated  from  the  2  GOM  tests.  In  this  (and  other  cases)  the  single  NS  test  suggests 
rather  lighter  damping;  indeed,  an  effectively  zero  value  of  pitch  damping  is  not  always  able  to  give 
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Model 

Description 

GOM 

N 

\S 

Cl 

Cs 

Cl 

Cs 

DF 

Base-case  model  with  diffraction  forces 

4.5 

1.6 

1.7 

.001 

DF/WDD 

DF  model  plus  wave-drift  damping 

3.3 

0.6 

.001 

.001 

DF/WDD/VF{u) 

DF/WDD  model  plus  viscous  forces  from 
undisturbed  waves 

4.0 

0.5 

2.5 

.001 

DF/WDD/VF(d) 

DF /WDD  model  plus  viscous  forces  from 
disturbed  waves 

6.5 

.001 

0.1 

.001 

Table  1:  Description  of  4  models,  and  consistent  damping  ratios  Ci  and  Cs  in  surge  and  pitch. 


as  narrow  a  spectral  bandwidth  as  observed.  Note,  however,  that  the  NS  test  includes  only  a  single 
hour,  making  narrow  bandwidths  more  difficult  to  estimate  than  over  the  combined,  two  hours  of 
GOM  tests. 

Our  second  model  (DF/WDD)  includes  both  diffraction  forces  and  wave-drift  damping.  This 
damping  force  is  proportional  both  to  the  structural  velocity  and  to  the  square  of  the  wave  am¬ 
plitude.  The  resulting  nonlinear  damping  will  tend  to  offset  large  slow  drift  forces,  and  perhaps 
reduce  the  extreme  peaks  of  the  surge  response.  We  may  expect  that  once  wave  drift  damping  is 
added,  we  require  lower  values  of  the  additional  modal  damping  Q  than  in  the  DF  model.  Table  1 
shows  that  this  is  indeed  the  case. 

Finally,  we  also  implement  two  models  that  include  viscous  drag  forces  as  well  as  diffraction  effects. 
These  differ  only  in  their  choice  of  wave  input:  one  uses  the  undisturbed  waves,  while  the  other 
uses  the  disturbed  waves  near  the  spar,  inferred  from  its  reported  heave  motions  and  the  air-gap 
(structure-to-wave  distance)  measurements.  Both  models  use  the  Morison’s  drag  term  with  C£>=0.6, 
and  Wheeler  stretching  (Wheeler,  1970)  to  integrate  effects  from  the  spar  bottom  to  the  free  surface. 
They  also  both  use  the  absolute  fluid  velocity;  relative  velocity  effects  are  assumed  reflected  through 
damping  terms. 

Numerical  Results 
Wave-Frequency  Response 

We  first  compare  the  wave-frequency  portions  of  the  predicted  and  observed  spar  responses.  These 
observed  portions  are  found  by  applying  a  high-pass  filter,  with  a  low-frequency  cutoff  of  .03  Hz., 
to  the  measured  displacement  histories.  Figure  5  shows  that  these  predictions  fairly  accurately 
predict  not  only  the  qualitative  response  behavior,  but  also  its  detailed  cycle-by-cycle  evolution  in 
all  three  tests.  This  suggests  that  our  models  accurately  reflect  first-order  wave  forces,  and  the  mass 
properties  of  the  spar.  (Because  slow-drift  forces  and  damping  do  not  affect  this  wave-frequency 
response,  all  four  of  our  models  predict  roughly  the  same  histories  in  Figure  5.  Thus,  for  clarity 
Figure  5  shows  predictions  only  for  the  simplest  (DF)  model.) 

Slow  Drift  Response 

Figure  6  compares  the  predicted  and  observed  moments  of  the  slow-drift  response.  Results  are 
shown  for  the  rms  values  (7i  and  as,  corresponding  to  frequency  ranges  (0-.006  Hz)  and  (.006-.03 
Hz),  and  for  the  total  mean  offset  (which  cannot  be  split  directly  into  surge  and  pitch  contributions). 
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Figure  6:  Predicted  response  mean  and  standard  deviations  vs  measured  results  in  the  three  tests. 
Standard  deviations  are  shown  for  both  surge  and  pitch  frequency  components. 


In  each  case  the  predicted  values  from  all  4  models  are  plotted  against  the  observed  value  from  that 
test.  Because  there  are  3  tests,  this  results  in  4  x  3=12  data  points  per  plot.  A  45-degree  line 
implying  perfect  agreement  is  also  shown. 

Focusing  first  on  the  diffraction-only  (DF)  model,  we  find  it  underpredicts  both  the  mean  and  rms 
in  all  but  one  of  the  9  observed  values.  (We  return  below  to  the  anomalous  case,  which  involves  the 
pitch  rms  <75  in  the  NS  test.)  Adding  wave  drift  damping  (DF/WDD)  will  not  change  the  mean, 
and  only  weakly  affects  the  rms  predictions  (again  excluding  the  NS  pitch  c£Lse).  This  suggests  the 
need  for  viscous  forces,  which  contribute  both  an  extra  mean  force  (due  to  the  asymmetric  effect 
of  wave  stretching)  and  a  slowly-varyiiig  drift  force.  As  Figure  6  shows,  the  VF(u)  and  VF(d) 
models — which  include  viscous  forces — generally  give  better  predictions  of  both  the  mean  and  rms 
levels.  Results  with  the  disturbed  wave  (VF(d))  generally  give  slightly  higher  responses — both  in 
mean  and  rms — than  those  using  the  undisturbed  wave  (VT'(u)).  Neither  the  VP(u)  nor  VF(d) 
model  seems  systematically  closer  to  the  observations;  however,  both  appear  superior  to  models 
that  exclude  viscous  forces  altogether. 

Returning  to  the  anamolous  pitch  response  in  the  NS  test,  note  from  Table  1  that  our  damping 
calibration  effectively  fails  in  this  case.  Although  each  of  the  4  models  was  assigned  only  minimal 
pitch  damping  (C5=.001),  all  of  these  predict  wider  spectral  bandwidths  than  that  observed  in  the 
NS  test.  Thus  the  predictive  models  here  are  not  “damping-tuned”  to  the  tests  as  in  the  other 
cases — and  the  pattern  of  the  4  model  predictions  for  as  in  the  NS  case  is  somewhat  arbitrary.  As 
to  why  the  bandwidth  mismatch  may  occur,  recall  the  increased  effect  of  limited  data  in  the  NS 
case:  here  the  bandwidth  estimation  uses  only  the  1  hour  test,  as  opposed  to  the  2  pooled  hours 
used  to  form  the  predictive  model  for  both  GOM  seastates. 

Total  Response  Histories 

Finally,  we  compare  the  observed  and  predicted  1-hour  histories  of  the  total  spar  displacement. 
Figures  7-9  show  these  histories  for  the  3  1-hour  tests.  All  figures  show  the  observed  displacement 
history  at  the  top,  while  3  of  the  4  corresponding  predictions  are  shown  beneath  (the  DF/WDD 


91 


model  is  omitted  for  clarity).  Recall  that  our  particular  interest  lies  in  predicting  not  only  overall 
response  statistics,  but  also  the  response  evolution  and  potential  mode-swapping  (e.g.,  Figure  4  for 
the  GOMl  case).  Figure  7  repeats  that  case,  and  it  is  notable  that  all  3  predicted  responses  show 
a  similar  trend,  toward  greater  surge  response,  in  the  second  half  of  the  GOMl  test.  Note  also 
that  for  extreme  response  events  (e.g.,  observed  response  above  20m),  all  of  the  predictions  show 
fairly  good  agreement.  The  greatest  deviations,  between  the  test  and  predictions,  seem  to  occur 
over  periods  of  relatively  low  response  amplitude  (e.g.,  times  t =0-1 000s,  1500-2300s). 

Figure  8  shows  similar  tendencies  for  the  GOM2  test.  Again  there  is  good  qualitative  agreement: 
both  the  tests  and  the  predictions  show  a  period  of  relatively  little  surge  (at  around  t=700-2000s), 
followed  by  a  marked  surge  increase  through  the  rest  of  the  hour.  Large  observed  responses  tend 
to  coincide  with  high  predicted  values.  The  magnitude  of  these  large  responses  is  not  as  well  pre¬ 
dicted,  however;  predictions  generally  underestimate  the  response  in  the  critical  high-surge  portion 
(i=2300-3600s).  In  contrast,  the  same  predictions  often  overestimate  response  in  the  earlier,  low- 
surge  segment  (^=700-2000s).  (This  potential  for  mode  swapping  and  compensating  errors  suggests 
the  need  here  to  compare  observations  and  predictions  through  time  histor)^  behavior,  not  merel}^ 
through  summar)^  statistics  in  the  time  or  frequency  domain.) 

Finally,  Figure  9  shows  NS  test  results.  Again  there  is  a  transition,  near  the  end  of  the  test,  which 
produces  the  largest  amplitude  responses  (indeed,  the  largest  offsets  among  the  3  tests).  It  is  caused 
here,  however,  by  an  increase  in  the  pitch  as  well  as  the  surge  component.  Note  that  despite  the 
potential  damping  mismatch  in  this  case,  the  models  follow  this  critical  portion  of  high-amplitude 
pitch  fairly  well  (from  about  t=2700s  on).  As  in  the  GOMl  case,  greater  deviations  between  models 
and  observations  occur  at  earlier  portions  of  the  historj^  involving  lower  amplitude  responses. 

Comparing  Responses  to  Simulated  and  Observed  Waves 

The  foregoing  results  show  how  well  various  models  can  predict  the  spar  response  in  the  OTRC 
tests,  based  on  the  corresponding  observed  wave  input.  Finally,  we  study  briefly  an  associated 
question:  are  the  observed  response  properties,  such  as  mode  swapping,  also  consistent  with  the 
response  of  the  spar  to  idealized,  Gaussian  simulations  of  (first-order)  random  waves?  If  not,  these 
observed  properties  may  perhaps  reflect  special  aspects  of  the  wave  input  in  the  wave  tank;  e.g., 
the  effect  of  its  finite  dimensions. 

As  earlier  noted,  both  surge  and  pitch  rms  components  of  the  observed  spar  motions  vary  notably 
over  periods  of  roughly  20-30  minutes.  We  quantify  this  rms  variation  b}^  (1)  splitting  the  response 
into  20-minute  segments;  (2)  calculating  the  rms  values,  ai...a„,  in  each  of  the  n  segments;  and  (3) 
forming  the  sample  mean  af/n  and  variance  “  ^)^/(^  !)•  We  focus  here  on  the 

two  GOM  tests,  yielding  2  hours  and  hence  n=6  20-minute  segments  with  associated  rms  values 
ai...aQ.  The  resulting  values  are  found  to  be 

Sff  =  1.22m  (surge) ;  Sa  =  0.50m  (pitch)  (3) 

For  comparison  we  simulate!  multiple  2-hour  spar  histories,  and  process  each  as  we  did  the  2-hour 
test  to  find  a  corresponding  s^r  estimate.  These  estimates  from  our  simulations  yield  the  following 

fThese  simulations  use  the  DF/WDD/VF(u)  model,  and  first-order  Gaussian  waves  are  simulated  from  a  JON- 
SVVAP  spectrum  with  i/,=14m,  Tp=14s,  and  7=2.  The  diffraction  analysis  internally  applies  second-order  corrections 
to  the  (assumed)  first-order  wave  input;  hence  the  Gaussian  model  is  consistent  here.  Drag  forces  for  this  model 
use  the  total  undisturbed  wave;  for  simplicity  we  use  the  Gaussian  waves  here  as  well.  Alternatively,  one  may  add 
second-order  wave  contributions  to  better  approximate  the  total  undisturbed  wave. 
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mean  £[5^7]  and  standard  deviation  D[5|y]: 

E[S(r]  =  L26m  (surge) ;  =  0.95m  (pitch)  (4) 

D[sa]  =  0.47m  (surge) ;  D[Sff]  =  0.45m  (pitch)  (5) 

Thus,  while  the  observed  modal  rms  values  may  seem  highly  variable,  our  simulations  show  similar 
or  still  prea^er  variability  (especially  in  pitch).  Note  also  that  the  D(.)  values  here  reflect  variability 
in  Sff  estimates  from  different  T=2  hour  segments.  (If  T  increases,  D{.)  should  decay  like  T"^/^.) 
Because  Eq.  3  uses  1  r=2  hour  segment,  these  D{.)  values  suggest  the  following  mean  ±l-sigma 
intervals  on  the  test  estimates:  s<y=L22±0.47  in  surge  and  Sey=0.50  ±  0.45  in  pitch.  As  even  these 
relatively  narrow,  1-sigma  confidence  intervals  include  the  average  simulation  results  (Eq.  4),  it  is 
difficult  to  find  statistically  significant  differences  between  the  tests  and  the  simulations. 

Conclusions 

Four  models  have  been  established  to  predict  the  along-wave  motions  of  a  spar  buoy  in  random  seas. 
These  have  been  implemented  and  compared  with  wave  tank  measurements  of  the  spar  displacement, 
at  a  reference  elevation  2:=54.8m  above  the  mean  water  level.  Results  are  shown  across  3  1-hour 
tests  of  100-year  extreme  wave  conditions.  Specific  methods  and  results  include  the  following: 

•  In  all  of  the  tests,  the  main  rms  contribution  comes  from  the  resonant  response  in  surge  and 
pitch  modes,  at  periods  of  roughly  330s  and  70s  respectively  (e.g..  Figure  3).  This  shows  the  need 
for  models  of  nonlinear  forces — diffraction,  drag  or  both — to  explain  not  only  the  mean  offset  but 
also  the  amplitude  of  slow-drift  oscillations  of  the  spar. 

•  The  significant  low-frequency  resonant  response  also  implies  the  need  for  accurate  estimates  of 
damping,  in  both  the  surge  and  pitch  modes  of  the  spar.  We  show  how  these  modal  dampings  can 
be  estimated  from  response  spectral  moments.  The  resulting  dampings  are  “consistent”  with  the 
other  features  of  the  model;  for  example,  the  explicit  addition  of  wave  drift  damping  (WDD)  is 
accompanied  by  lower  levels  of  the  remaining  damping  in  the  model  (Table  1). 

•  The  wave-frequency  response  has  been  found  to  be  fairly  well  predicted  across  all  3  tests  (Fig¬ 
ure  5).  This  reflects  the  modelling  adequacy  of  linear  diffraction  forces  and  the  spar’s  mass  proper¬ 
ties.  Regarding  slow-drift  response,  models  that  include  only  diffraction  forces  generally  underesti¬ 
mate  both  the  mean  and  rms  response  levels  (Figure  6).  To  address  this,  we  introduce  additional 
models  that  include  viscous  forces,  based  on  either  the  undisturbed  (far-field)  wave  or  the  actual 
disturbed  wave  in  the  presence  of  the  spar.  While  it  is  difficult  to  conclude  which  of  these  is  generally 
more  accurate,  both  appear  superior  to  models  that  exclude  viscous  forces  altogether. 

•  The  observed  responses  display  considerable  “mode-swapping”  between  surge  and  pitch  modes 
(e.g.,  Figure  4).  Figures  7-9  show  that  our  predictive  models,  which  use  the  observed  wave  and 
its  underlying  first-order  components,  can  produce  qualitatively  similar  behavior.  They  generally 
follow  the  observed  trend,  in  all  three  tests,  toward  larger  amplitude  responses  near  the  end  of  the 
hour.  This  trend  manifests  itself  in  the  Gulf  of  Mexico  tests  by  a  late  increase  in  surge-induced 
response  (Figures  7-8),  and  in  the  North  Sea  test  by  enhanced  pitch  response  as  well  (Figure  9). 

•  While  the  modal  rms  values  in  the  tests  appear  rather  variable,  long  simulations  with  Gaussian 
(first-order)  waves  show  similar  or  still  greater  variability  (Eqs.  3-5).  From  the  limited  2-hour 
duration  of  GOM  tests,  it  is  difficult  to  find  statistically  significant  differences  between  these  tests 
and  the  simulations. 
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Figure  9:  Combined  (total)  surge  response  time  history  for  NS:  prediction  vs.  measurement. 
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